sentieon使用

  • gvcf_joincalling

vi test.sh
ls *.gz >input.setieon.list
export SENTIEON_LICENSE=192.168.102.80:8990
cat input.setieon.list | /home/yzhou/soft/sentieon-genomics-202112.07/bin/sentieon driver -r /home/yzhou/refs/pig/Sus_scrofa.Sscrofa11.1-106/Sus_scrofa.Sscrofa11.1.dna.toplevel.fa -t 20 --algo GVCFtyper -d /home/yzhou/refs/pig/Sus_scrofa.Sscrofa11.1-106/sus_scrofa.vcf.gz /data1/yzhou/ZZQ2003/typedvcf/combZZQ2003_chr13_typed.vcf.gz -
  • 单个个体SNPcalling

cd /data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata
mkdir DNAseq
vi  test01.sh
#!/bin/sh
# Copyright (c) 2016-2020 Sentieon Inc. All rights reserved
# *******************************************
# Script to perform DNA seq variant calling
# using a single sample with fastq files
# named 1.fastq.gz and 2.fastq.gz
# *******************************************
set -eu
# Update with the fullpath location of your sample fastq
SM="test01" #sample name
RGID="rg_$SM" #read group ID
PL="ILLUMINA" #or other sequencing platform
FASTQ_FOLDER="/data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata"
FASTQ_1="$FASTQ_FOLDER/${SM}_1_val_1.fq.gz"
FASTQ_2="$FASTQ_FOLDER/${SM}_2_val_2.fq.gz"  #If using Illumina paired data
# Update with the location of the reference data files
FASTA_DIR="/home/jychu/home/refs/pig/Sus_scrofa.Sscrofa11.1-106"
FASTA="$FASTA_DIR/Sus_scrofa.Sscrofa11.1.dna.toplevel.fa"
KNOWN_DBSNP="$FASTA_DIR/sus_scrofa.vcf.gz"
# Update with the location of the Sentieon software package and license file
SENTIEON_INSTALL_DIR=/home/yzhou/soft/sentieon-genomics-202112.07
export SENTIEON_LICENSE=192.168.102.80:8990  #or using licsrvr: c1n11.sentieon.com:5443
# Other settings
NT=1  #number of threads to use in computation, set to number of cores in the server
START_DIR="$FASTQ_FOLDER/DNAseq" #Determine where the output files will be stored
# You do not need to modify any of the lines below unless you want to tweak the pipeline
# ************************************************************************************************************************************************************************
# ******************************************
# 0. Setup
# ******************************************
WORKDIR="$START_DIR/${SM}" 
mkdir -p $WORKDIR
LOGFILE=$WORKDIR/run.log
exec >$LOGFILE 2>&1
cd $WORKDIR
# ******************************************
# 1. Mapping reads with BWA-MEM, sorting
# ******************************************
#The results of this call are dependent on the number of threads used. To have number of threads independent results, add chunk size option -K 10000000 
( $SENTIEON_INSTALL_DIR/bin/sentieon bwa mem -R "@RG\tID:$RGID\tSM:$SM\tPL:$PL" \
    -t $NT -K 10000000 $FASTA $FASTQ_1 $FASTQ_2 || { echo -n 'BWA error'; exit 1; } ) | \
    $SENTIEON_INSTALL_DIR/bin/sentieon util sort -r $FASTA -o sorted.bam -t $NT \
    --sam2bam -i - || { echo "Alignment failed"; exit 1; }
# ******************************************
# 2. Metrics
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i sorted.bam \
    --algo MeanQualityByCycle mq_metrics.txt --algo QualDistribution qd_metrics.txt \
    --algo GCBias --summary gc_summary.txt gc_metrics.txt --algo AlignmentStat \
    --adapter_seq '' aln_metrics.txt --algo InsertSizeMetricAlgo is_metrics.txt || \
    { echo "Metrics failed"; exit 1; }

$SENTIEON_INSTALL_DIR/bin/sentieon plot GCBias -o gc-report.pdf gc_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot QualDistribution -o qd-report.pdf qd_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot MeanQualityByCycle -o mq-report.pdf mq_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot InsertSizeMetricAlgo -o is-report.pdf is_metrics.txt

# ******************************************
# 3. Remove Duplicate Reads. It is possible
# to remove instead of mark duplicates
# by adding the --rmdup option in Dedup
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo LocusCollector \
    --fun score_info score.txt || { echo "LocusCollector failed"; exit 1; }
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo Dedup \
    --score_info score.txt --metrics dedup_metrics.txt deduped.bam || \
    { echo "Dedup failed"; exit 1; }
# ******************************************
# 5. Base recalibration
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam --algo QualCal \
    -k $KNOWN_DBSNP recal_data.table
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
    --algo QualCal -k $KNOWN_DBSNP  recal_data.table.post
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT --algo QualCal --plot \
    --before recal_data.table --after recal_data.table.post recal.csv
$SENTIEON_INSTALL_DIR/bin/sentieon plot QualCal -o recal_plots.pdf recal.csv
# ******************************************
# 6b. HC Variant caller
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
    --algo Haplotyper -d $KNOWN_DBSNP ${SM}-hc.vcf.gz || \
    { echo "Haplotyper failed"; exit 1; }
# ******************************************
# 5b. ReadWriter to output recalibrated bam
# This stage is optional as variant callers
# can perform the recalibration on the fly
# using the before recalibration bam plus
# the recalibration table
# ******************************************
# $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table  --algo ReadWriter recaled.bam || { echo "ReadWriter failed"; exit 1; }
for i in `ls *_1_*|tr "_" "\t"|cut -f1`;do sed "s/test01/$i/g"  test01.sh >SNP_CALL$i.sh;done
for i in `ls *_1_*|tr "_" "\t"|cut -f1`;do  chmod a+x SNP_CALL$i.sh;done
  • 大群体joincalling

cd /data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata
mkdir DNAseq
vi  test01.sh
#!/bin/sh
# Copyright (c) 2016-2020 Sentieon Inc. All rights reserved
# *******************************************
# Script to perform DNA seq variant calling
# using a single sample with fastq files
# named 1.fastq.gz and 2.fastq.gz
# *******************************************
cd  /data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata
set -eu
# Update with the fullpath location of your sample fastq
SM="test01" #sample name
RGID="rg_$SM" #read group ID
PL="ILLUMINA" #or other sequencing platform
FASTQ_FOLDER="/data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata"
FASTQ_1="$FASTQ_FOLDER/${SM}_1_val_1.fq.gz"
FASTQ_2="$FASTQ_FOLDER/${SM}_2_val_2.fq.gz"  #If using Illumina paired data
# Update with the location of the reference data files
FASTA_DIR="/home/jychu/home/refs/pig/Sus_scrofa.Sscrofa11.1-106"
FASTA="$FASTA_DIR/Sus_scrofa.Sscrofa11.1.dna.toplevel.fa"
KNOWN_DBSNP="$FASTA_DIR/sus_scrofa.vcf.gz"
# Update with the location of the Sentieon software package and license file
SENTIEON_INSTALL_DIR=/home/yzhou/soft/sentieon-genomics-202112.07
export SENTIEON_LICENSE=192.168.102.80:8990  #or using licsrvr: c1n11.sentieon.com:5443
# Other settings
NT=1  #number of threads to use in computation, set to number of cores in the server
START_DIR="$FASTQ_FOLDER/DNAseq" #Determine where the output files will be stored
# You do not need to modify any of the lines below unless you want to tweak the pipeline
# ************************************************************************************************************************************************************************
# ******************************************
# 0. Setup
# ******************************************
WORKDIR="$START_DIR/${SM}" 
mkdir -p $WORKDIR
LOGFILE=$WORKDIR/run.log
exec >$LOGFILE 2>&1
cd $WORKDIR
# ******************************************
# 1. Mapping reads with BWA-MEM, sorting
# ******************************************
#The results of this call are dependent on the number of threads used. To have number of threads independent results, add chunk size option -K 10000000 
( $SENTIEON_INSTALL_DIR/bin/sentieon bwa mem -R "@RG\tID:$RGID\tSM:$SM\tPL:$PL" \
    -t $NT -K 10000000 $FASTA $FASTQ_1 $FASTQ_2 || { echo -n 'BWA error'; exit 1; } ) | \
    $SENTIEON_INSTALL_DIR/bin/sentieon util sort -r $FASTA -o sorted.bam -t $NT \
    --sam2bam -i - || { echo "Alignment failed"; exit 1; }
# ******************************************
# 2. Metrics
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i sorted.bam \
    --algo MeanQualityByCycle mq_metrics.txt --algo QualDistribution qd_metrics.txt \
    --algo GCBias --summary gc_summary.txt gc_metrics.txt --algo AlignmentStat \
    --adapter_seq '' aln_metrics.txt --algo InsertSizeMetricAlgo is_metrics.txt || \
    { echo "Metrics failed"; exit 1; }

$SENTIEON_INSTALL_DIR/bin/sentieon plot GCBias -o gc-report.pdf gc_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot QualDistribution -o qd-report.pdf qd_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot MeanQualityByCycle -o mq-report.pdf mq_metrics.txt
$SENTIEON_INSTALL_DIR/bin/sentieon plot InsertSizeMetricAlgo -o is-report.pdf is_metrics.txt

# ******************************************
# 3. Remove Duplicate Reads. It is possible
# to remove instead of mark duplicates
# by adding the --rmdup option in Dedup
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo LocusCollector \
    --fun score_info score.txt || { echo "LocusCollector failed"; exit 1; }
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo Dedup \
    --score_info score.txt --metrics dedup_metrics.txt deduped.bam || \
    { echo "Dedup failed"; exit 1; }
# ******************************************
# 5. Base recalibration
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam --algo QualCal \
    -k $KNOWN_DBSNP recal_data.table
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
    --algo QualCal -k $KNOWN_DBSNP  recal_data.table.post
$SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT --algo QualCal --plot \
    --before recal_data.table --after recal_data.table.post recal.csv
$SENTIEON_INSTALL_DIR/bin/sentieon plot QualCal -o recal_plots.pdf recal.csv
 # ******************************************
 # 6b. HC Variant caller for GVCF
 # ******************************************
 $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
      --algo Haplotyper -d $KNOWN_DBSNP --emit_conf=30 --call_conf=30 --emit_mode gvcf \
      ${SM}-hc.g.vcf.gz || { echo "Haplotyper failed"; exit 1; }

for i in `ls *_1_*|tr "_" "\t"|cut -f1`;do sed "s/test01/$i/g"  test01.sh >SNP_CALL$i.sh;done
for i in `ls *_1_*|tr "_" "\t"|cut -f1`;do  chmod a+x SNP_CALL$i.sh;done
  • 小群体joincalling

vi joincalling.sh
 #!/bin/sh
# Copyright (c) 2016-2020 Sentieon Inc. All rights reserved
# *******************************************
# Script to perform joint calling on 3 samples
# with fastq files named sample<i>_1.fastq.gz
# and sample<i>_2.fastq.gz
# *******************************************
set -eu
# Update with the fullpath location of your sample fastq
SM_LIST=$(ls *_1_val_1.fq.gz|tr "_" "\t"|cut -f1) # list of sample names
PL="ILLUMINA" #or other sequencing platform
FASTQ_FOLDER="/data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata"
FASTQ_1_SUFFIX="1_val_1.fq.gz"
FASTQ_2_SUFFIX="2_val_2.fq.gz"
# Update with the location of the reference data files
FASTA_DIR="/home/jychu/home/refs/pig/Sus_scrofa.Sscrofa11.1-106"
FASTA="$FASTA_DIR/Sus_scrofa.Sscrofa11.1.dna.toplevel.fa"
KNOWN_DBSNP="$FASTA_DIR/sus_scrofa.vcf.gz"
# Update with the location of the Sentieon software package and license file
SENTIEON_INSTALL_DIR=/home/yzhou/soft/sentieon-genomics-202112.07
export SENTIEON_LICENSE=192.168.102.80:8990  #or using licsrvr: c1n11.sentieon.com:5443
# Other settings
NT=40  #number of threads to use in computation, set to number of cores in the server
START_DIR="/data1/data_public/WGS/pig/XinJiang382/2.5G/rawdata/cleandata/DNAseq_joint" #Determine where the output files will be stored
# You do not need to modify any of the lines below unless you want to tweak the pipeline
# ************************************************************************************************************************************************************************
# ******************************************
# 0. Setup
# ******************************************
WORKDIR="$START_DIR/joint_call"
mkdir -p $WORKDIR
LOGFILE=$WORKDIR/run.log
exec >$LOGFILE 2>&1
# ******************************************
# 0. Process all samples independently
# ******************************************
GVCF_INPUTS=""
for SM in $SM_LIST; do
  RGID="rg_$SM"
  mkdir $WORKDIR/$SM
  cd $WORKDIR/$SM
  # ******************************************
  # 1. Mapping reads with BWA-MEM, sorting
  # ******************************************
  ( $SENTIEON_INSTALL_DIR/bin/sentieon bwa mem -R "@RG\tID:$RGID\tSM:$SM\tPL:$PL" \
      -t $NT -K 10000000 $FASTA $FASTQ_FOLDER/${SM}_$FASTQ_1_SUFFIX \
      $FASTQ_FOLDER/${SM}_$FASTQ_2_SUFFIX || { echo -n 'bwa error'; exit 1; } ) | \
      $SENTIEON_INSTALL_DIR/bin/sentieon util sort -r $FASTA -o sorted.bam -t $NT --sam2bam -i - || \
      { echo "Alignment failed"; exit 1; }
  # ******************************************
  # 2. Metrics
  # ******************************************
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i sorted.bam \
      --algo MeanQualityByCycle mq_metrics.txt --algo QualDistribution qd_metrics.txt \
      --algo GCBias --summary gc_summary.txt gc_metrics.txt --algo AlignmentStat \
      --adapter_seq '' aln_metrics.txt --algo InsertSizeMetricAlgo is_metrics.txt || \
      { echo "Metrics failed"; exit 1; }
  $SENTIEON_INSTALL_DIR/bin/sentieon plot GCBias -o gc-report.pdf gc_metrics.txt
  $SENTIEON_INSTALL_DIR/bin/sentieon plot QualDistribution -o qd-report.pdf qd_metrics.txt
  $SENTIEON_INSTALL_DIR/bin/sentieon plot MeanQualityByCycle -o mq-report.pdf mq_metrics.txt
  $SENTIEON_INSTALL_DIR/bin/sentieon plot InsertSizeMetricAlgo -o is-report.pdf is_metrics.txt
 
  # ******************************************
  # 3. Remove Duplicate Reads. It is possible
  # to remove instead of mark duplicates
  # by adding the --rmdup option in Dedup
  # ******************************************
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo LocusCollector \
      --fun score_info score.txt || { echo "LocusCollector failed"; exit 1; }

  $SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT -i sorted.bam --algo Dedup \
      --score_info score.txt --metrics dedup_metrics.txt deduped.bam || \
      { echo "Dedup failed"; exit 1; }
  # ******************************************
  # 2a. Coverage metrics
  # ******************************************
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam \
      --algo CoverageMetrics coverage_metrics || { echo "CoverageMetrics failed"; exit 1; }
  # ******************************************
  # 5. Base recalibration
  # ******************************************
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam --algo QualCal \
      -k $KNOWN_DBSNP  recal_data.table
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
      --algo QualCal -k $KNOWN_DBSNP  recal_data.table.post
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -t $NT --algo QualCal --plot \
      --before recal_data.table --after recal_data.table.post recal.csv
  $SENTIEON_INSTALL_DIR/bin/sentieon plot QualCal -o recal_plots.pdf recal.csv
  # ******************************************
  # 6b. HC Variant caller for GVCF
  # ******************************************
  $SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT -i deduped.bam -q recal_data.table \
      --algo Haplotyper -d $KNOWN_DBSNP --emit_conf=30 --call_conf=30 --emit_mode gvcf \
      output-hc.g.vcf.gz || { echo "Haplotyper failed"; exit 1; }
  GVCF_INPUTS="$GVCF_INPUTS -v $WORKDIR/$SM/output-hc.g.vcf.gz"
done
# ******************************************
# Perform the joint calling 
# ******************************************
$SENTIEON_INSTALL_DIR/bin/sentieon driver -r $FASTA -t $NT --algo GVCFtyper $GVCF_INPUTS \
    -d $KNOWN_DBSNP $WORKDIR/output-jc.vcf.gz || { echo "GVCFtyper failed"; exit 1; }
chmod a+x joincalling.sh
nohup ./joincalling.sh &
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 211,884评论 6 492
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,347评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,435评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,509评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,611评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,837评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,987评论 3 408
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,730评论 0 267
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,194评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,525评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,664评论 1 340
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,334评论 4 330
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,944评论 3 313
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,764评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,997评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,389评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,554评论 2 349

推荐阅读更多精彩内容