VariantBam
VariantBam copied to clipboard
walaj
Soft-clipping with linked-region filter
Hi, first -- thank you for your software! I am working with paired end sequencing on a pool of fragments that contain a large proportion of fragments shorter than the reads. This leads to soft-clipping in the (proper) alignment. Unfortunately, when filtering paired reads using the -l (linked region) flag, proper read pairs get included that do not overlap the given locus.
I've provided a MWE:
## github.sam ##
@HD VN:1.5 SO:coordinate
@SQ SN:NC_000001.11 LN:248956422
@RG ID:S3_L002 SM:S3 PL:Illumina CN:UNC
@RG ID:S3_L003 SM:S3 PL:Illumina CN:UNC
@PG ID:bwa PN:bwa VN:0.7.15-r1140 CL:bwa mem -t 8 -R @RG\tID:S3_L002\tSM:S3\tPL:Illumina\tCN:UNC /projects/sequence_analysis/vol4/bpow/references/GCF_000001405.26_GRCh38_genomic.fna inputs/S3/L002/S3_L002_R2.fastq.gz inputs/S3/L002/S3_L002_R1.fastq.gz
@PG ID:MarkDuplicates VN:2.9.2-SNAPSHOT CL:picard.sam.markduplicates.MarkDuplicates INPUT=[mapped/S3.L002.bam] OUTPUT=markdup/S3.L002.bam METRICS_FILE=markdup/S3.L002.metrics.txt REMOVE_DUPLICATES=false TMP_DIR=[/scratch/dfiler] CREATE_INDEX=true MAX_SEQUENCES_FOR_DISK_READ_ENDS_MAP=50000 MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 SORTING_COLLECTION_SIZE_RATIO=0.25 TAG_DUPLICATE_SET_MEMBERS=false REMOVE_SEQUENCING_DUPLICATES=false TAGGING_POLICY=DontTag ASSUME_SORTED=false DUPLICATE_SCORING_STRATEGY=SUM_OF_BASE_QUALITIES PROGRAM_RECORD_ID=MarkDuplicates PROGRAM_GROUP_NAME=MarkDuplicates READ_NAME_REGEX=<optimized capture of last three ':' separated fields as numeric values> OPTICAL_DUPLICATE_PIXEL_DISTANCE=100 VERBOSITY=INFO QUIET=false VALIDATION_STRINGENCY=STRICT COMPRESSION_LEVEL=5 MAX_RECORDS_IN_RAM=500000 CREATE_MD5_FILE=false GA4GH_CLIENT_SECRETS=client_secrets.json PN:MarkDuplicates
@PG ID:bwa-4A25094 PN:bwa VN:0.7.15-r1140 CL:bwa mem -t 8 -R @RG\tID:S3_L003\tSM:S3\tPL:Illumina\tCN:UNC /projects/sequence_analysis/vol4/bpow/references/GCF_000001405.26_GRCh38_genomic.fna inputs/S3/L003/S3_L003_R1.fastq.gz inputs/S3/L003/S3_L003_R2.fastq.gz
@PG ID:MarkDuplicates-4B2F2425 VN:2.9.2-SNAPSHOT CL:picard.sam.markduplicates.MarkDuplicates INPUT=[mapped/S3.L003.bam] OUTPUT=markdup/S3.L003.bam METRICS_FILE=markdup/S3.L003.metrics.txt REMOVE_DUPLICATES=false TMP_DIR=[/scratch/dfiler] CREATE_INDEX=true MAX_SEQUENCES_FOR_DISK_READ_ENDS_MAP=50000 MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 SORTING_COLLECTION_SIZE_RATIO=0.25 TAG_DUPLICATE_SET_MEMBERS=false REMOVE_SEQUENCING_DUPLICATES=false TAGGING_POLICY=DontTag ASSUME_SORTED=false DUPLICATE_SCORING_STRATEGY=SUM_OF_BASE_QUALITIES PROGRAM_RECORD_ID=MarkDuplicates PROGRAM_GROUP_NAME=MarkDuplicates READ_NAME_REGEX=<optimized capture of last three ':' separated fields as numeric values> OPTICAL_DUPLICATE_PIXEL_DISTANCE=100 VERBOSITY=INFO QUIET=false VALIDATION_STRINGENCY=STRICT COMPRESSION_LEVEL=5 MAX_RECORDS_IN_RAM=500000 CREATE_MD5_FILE=false GA4GH_CLIENT_SECRETS=client_secrets.json PN:MarkDuplicates
K00269:198:H55LWBBXY:2:2222:28392:34618 83 NC_000001.11 1180725 60 29S121M = 1180727 -119 ACTCTTTCCCTACACGACGCTCTTCCGATCTCCAGGAGAGCAGCTGCTGTACCAGCTTCCCAACAACAAGCTCCTCACCACCAAGATCGGGCTGCTCAGCACCCTTCGGGGACGGGCACGGGCCATGAGCAAGGCCAGCAAGGTGCCGGG -JJJJJFJJJJJJJJJFJJJJJJJJJJJJJFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFFFAA MD:Z:121 NM:i:0 AS:i:121 XS:i:0 RG:Z:S3_L002 PG:Z:MarkDuplicates
K00269:198:H55LWBBXY:2:2222:28392:34618 163 NC_000001.11 1180727 60 119M31S = 1180725 119 CCAGGAGAGCAGCTGCTGTACCAGCTTCCCAACAACAAGCTCCTCACCACCAAGATCGGGCTGCTCAGCACCCTTCGGGGACGGGCACGGGCCATGAGCAAGGCCAGCAAGGTGCCGGGAGATCGGAAGAGCACACGTCTGAACTCCAGT AAFFFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFJJJJJJ MD:Z:119 NM:i:0 AS:i:119 XS:i:19 RG:Z:S3_L002 PG:Z:MarkDuplicates
## github.vcf ##
##fileformat=VCFv4.2
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT SMP1 SMP2
NC_000001.11 1180851 NC_000001.11:1180851_T/C T C 20450.4 . AB=0.471229;ABP=12.2568;AC=1;AF=0.25;AN=4;AO=831;CIGAR=9M1X;DP=133;DPB=8557.7;DPRA=1.09321;EPP=98.3382;EPPR=1443.91;GTI=0;LEN=1;MEANALT=2.60976;MQM=60;MQMR=59.9988;NS=126;NUMALT=3;ODDS=13.1657;PAIRED=0.998797;PAIREDR=0.997442;PAO=430.75;PQA=14556.7;PQR=14751.7;PRO=442.083;QA=29918;QR=245335;RO=6646;RPL=612;RPP=406.598;RPPR=3656.66;RPR=219;RUN=1;SAF=117;SAP=934.337;SAR=714;SRF=742;SRP=8709.24;SRR=5904;TYPE=snp;technology.Illumina=1 GT:AD:AO:DP:GQ:PL:QA:QR:RO 0/0:76,0:0:78:99:0,436,2669:0:2983:76 0/1:26,27:27:55:99:594,0,542:1065:1007:26
$ variant github.sam -l github.vcf
@HD VN:1.5 SO:coordinate
@SQ SN:NC_000001.11 LN:248956422
@RG ID:S3_L002 SM:S3 PL:Illumina CN:UNC
@RG ID:S3_L003 SM:S3 PL:Illumina CN:UNC
@PG ID:bwa PN:bwa VN:0.7.15-r1140 CL:bwa mem -t 8 -R @RG\tID:S3_L002\tSM:S3\tPL:Illumina\tCN:UNC /projects/sequence_analysis/vol4/bpow/references/GCF_000001405.26_GRCh38_genomic.fna inputs/S3/L002/S3_L002_R2.fastq.gz inputs/S3/L002/S3_L002_R1.fastq.gz
@PG ID:MarkDuplicates VN:2.9.2-SNAPSHOT CL:picard.sam.markduplicates.MarkDuplicates INPUT=[mapped/S3.L002.bam] OUTPUT=markdup/S3.L002.bam METRICS_FILE=markdup/S3.L002.metrics.txt REMOVE_DUPLICATES=false TMP_DIR=[/scratch/dfiler] CREATE_INDEX=true MAX_SEQUENCES_FOR_DISK_READ_ENDS_MAP=50000 MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 SORTING_COLLECTION_SIZE_RATIO=0.25 TAG_DUPLICATE_SET_MEMBERS=false REMOVE_SEQUENCING_DUPLICATES=false TAGGING_POLICY=DontTag ASSUME_SORTED=false DUPLICATE_SCORING_STRATEGY=SUM_OF_BASE_QUALITIES PROGRAM_RECORD_ID=MarkDuplicates PROGRAM_GROUP_NAME=MarkDuplicates READ_NAME_REGEX=<optimized capture of last three ':' separated fields as numeric values> OPTICAL_DUPLICATE_PIXEL_DISTANCE=100 VERBOSITY=INFO QUIET=false VALIDATION_STRINGENCY=STRICT COMPRESSION_LEVEL=5 MAX_RECORDS_IN_RAM=500000 CREATE_MD5_FILE=false GA4GH_CLIENT_SECRETS=client_secrets.json PN:MarkDuplicates
@PG ID:bwa-4A25094 PN:bwa VN:0.7.15-r1140 CL:bwa mem -t 8 -R @RG\tID:S3_L003\tSM:S3\tPL:Illumina\tCN:UNC /projects/sequence_analysis/vol4/bpow/references/GCF_000001405.26_GRCh38_genomic.fna inputs/S3/L003/S3_L003_R1.fastq.gz inputs/S3/L003/S3_L003_R2.fastq.gz
@PG ID:MarkDuplicates-4B2F2425 VN:2.9.2-SNAPSHOT CL:picard.sam.markduplicates.MarkDuplicates INPUT=[mapped/S3.L003.bam] OUTPUT=markdup/S3.L003.bam METRICS_FILE=markdup/S3.L003.metrics.txt REMOVE_DUPLICATES=false TMP_DIR=[/scratch/dfiler] CREATE_INDEX=true MAX_SEQUENCES_FOR_DISK_READ_ENDS_MAP=50000 MAX_FILE_HANDLES_FOR_READ_ENDS_MAP=8000 SORTING_COLLECTION_SIZE_RATIO=0.25 TAG_DUPLICATE_SET_MEMBERS=false REMOVE_SEQUENCING_DUPLICATES=false TAGGING_POLICY=DontTag ASSUME_SORTED=false DUPLICATE_SCORING_STRATEGY=SUM_OF_BASE_QUALITIES PROGRAM_RECORD_ID=MarkDuplicates PROGRAM_GROUP_NAME=MarkDuplicates READ_NAME_REGEX=<optimized capture of last three ':' separated fields as numeric values> OPTICAL_DUPLICATE_PIXEL_DISTANCE=100 VERBOSITY=INFO QUIET=false VALIDATION_STRINGENCY=STRICT COMPRESSION_LEVEL=5 MAX_RECORDS_IN_RAM=500000 CREATE_MD5_FILE=false GA4GH_CLIENT_SECRETS=client_secrets.json PN:MarkDuplicates
K00269:198:H55LWBBXY:2:2222:28392:34618 83 NC_000001.11 1180725 60 29S121M = 1180727 -119 ACTCTTTCCCTACACGACGCTCTTCCGATCTCCAGGAGAGCAGCTGCTGTACCAGCTTCCCAACAACAAGCTCCTCACCACCAAGATCGGGCTGCTCAGCACCCTTCGGGGACGGGCACGGGCCATGAGCAAGGCCAGCAAGGTGCCGGG -JJJJJFJJJJJJJJJFJJJJJJJJJJJJJFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFFFAA MD:Z:121 NM:i:0 AS:i:121 XS:i:0 RG:Z:S3_L002 PG:Z:MarkDuplicates
K00269:198:H55LWBBXY:2:2222:28392:34618 163 NC_000001.11 1180727 60 119M31S = 1180725 119 CCAGGAGAGCAGCTGCTGTACCAGCTTCCCAACAACAAGCTCCTCACCACCAAGATCGGGCTGCTCAGCACCCTTCGGGGACGGGCACGGGCCATGAGCAAGGCCAGCAAGGTGCCGGGAGATCGGAAGAGCACACGTCTGAACTCCAGT AAFFFJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJJFJJJJJJ MD:Z:119 NM:i:0 AS:i:119 XS:i:19 RG:Z:S3_L002 PG:Z:MarkDuplicates
For illustration, using ASCIIGenome:
$ ASCIIGenome github-out.sam github.vcf
github-out.sam@2; Reads: 2
ctccaggagagcagctgctgtaccagcttcccaacaacaagctcctcaccaccaagatcgggctgctcagcacccttcggggacgggcacgggccatgagcaaggccagcaaggtgccggg
CCAGGAGAGCAGCTGCTGTACCAGCTTCCCAACAACAAGCTCCTCACCACCAAGATCGGGCTGCTCAGCACCCTTCGGGGACGGGCACGGGCCATGAGCAAGGCCAGCAAGGTGCCGGG
github.vcf#3; N: 1 C
SMP1 .
SMP2 E
1180725 1180735 1180745 1180755 1180765 1180775 1180785 1180795 1180805 1180815 1180825 1180835 1180845 11
0 .08 .15 .23 .31 .38 .46 .53 .61 .69 .76 .84 .92 .9
NC_000001.11:1180725-1180856; 132 bp; 1.0 bp/char /\
You can see the mapped portion of the reads do not cover the locus of interest.
