29 September, 2020
SNPsplit is an allele-specific alignment sorter which is designed to read alignment files in SAM/ BAM format and determine the allelic origin of reads that cover known SNP positions. For this to work a library must have been aligned to a genome which had all SNP positions masked by the ambiguity base ’N’, and aligned using aligners that are capable of using a reference genome which contains ambiguous nucleobases. Examples of supported alignment programs are Bowtie2, HISAT2, Bismark, HiCUP, TopHat or STAR (for some tips using HISAT2 or STAR alignments please see below). In addition, a list of all known SNP positions between the two different genomes must be provided using the option --snp_file. SNP information to generate N-masked genomes needs to be acquired elsewhere, e.g. for different strains of mice you can find variant call files at the Mouse Genomes Project page at http://www.sanger.ac.uk/resources/mouse/genomes/.
SNPsplit offers a separate genome preparation step that allows you to generate N-masked (or fully incorporated) SNP genomes for single or dual hybrid strains for all strains of the Mouse Genomes Project. Please see below for further details.
SNPsplit operates in two stages:
I) SNPsplit-tag: SNPsplit analyses reads (single-end mode) or read pairs (paired-end mode) for overlaps with known SNP positions, and writes out a tagged BAM file in the same order as the original file. Unsorted paired-end files are sorted by name first.
II) SNPsplit-sort: The tagged BAM file is read in again and sorted into allele-specific files. This process may also be run as a stand-alone module on tagged BAM files (tag2sort).
The SNPsplit-tag module determines whether a read can be assigned to a certain allele and appends an additional optional field XX:Z: to each read. The tag can be one of the following:
XX:Z:UA - Unassigned
XX:Z:G1 - Genome 1-specific
XX:Z:G2 - Genome 2-specific
XX:Z:CF - Conflicting
The SNPsplit-sort module tag2sort reads in the tagged BAM file and sorts the reads (or read pairs) according to their XX:Z: tag (or the combination of tags for paired-end or Hi-C reads) into sub-files.
Optional. If the supplied file is a SAM file it will first be converted to BAM format (using samtools view).
Paired-end files might require the input file to be sorted by read ID before continuing with the allele-tagging (Read 1 and Read 2 of a pair are expected to follow each other in the input BAM file). Unless specifically stated, paired-end BAM files will be sorted by read name (using samtools sort -n; output file ending in .sortedByName.bam). For files that already contain R1 and R2 on two consecutive lanes, the sorting step may be skipped using the option --no_sort. Single-end files or Hi-C files generated by HiCUP do not require sorting.
SNP positions are read in from the SNP file (which may be GZIP compressed (ending in .gz) or plain text files). The SNP file is expected to be in the following format (tab-delimited):
ID Chr Position SNP value Ref/SNP
18819008 5 48794752 1 C/T
40491905 11 63643453 1 A/G
44326884 12 96627819 1 T/A
Only the information contained in fields 'Chr (Chromosome)', 'Position' and 'Ref/SNP' base are being used for analysis. The genome containing the 'Ref' base is used for 'genome 1 specific reads (G1)', the genome containing the 'SNP' base for 'genome 2 specific reads (G2)'. If reads do not overlap any SNP positions they are considered 'Unassigned (UA)', i.e. they are not informative for one allele or another. In the rare case that a read contains both genome 1- and genome 2-specific base(s), or that the SNP position was deleted the read is regarded as 'Conflicting (CF)'.
It is probably noteworthy that the determination of overlaps correctly handles the CIGAR operations M (match), D (deletion in the read), I (insertion in the read) and N (skipped regions, used for splice mapping by e.g. TopHat or HISAT2). Other CIGAR operations are currently not supported. This is the reason why SNPsplit does not support soft-clipping.
Upon completion, a small allele-specific tagging report is printed to screen and to a report file (.SNPsplit_report.txt) for archiving purposes.
Once the tagging has completed, the tag2sort module reads in the tagged BAM file and sorts it into various sub-files according to their XX:Z: tag. Both single and paired-end files are sorted into the following four categories:
UA - Unassigned
G1 - Genome 1-specific
G2 - Genome 2-specific
CF - Conflicting
Files with conflicting SNP information (tag XX:Z:CF) are not written out by default.
Upon completion, an allele-specific sorting report is printed out on screen and to a report file for archiving purposes (*.SNPsplit_sort.txt). If the sorting was launched by SNPsplit and not run stand-alone (as tag2sort) the sorting report will also be written into the main SNPsplit report (*.SNPsplit_report.txt).
In paired-end mode, both reads are used for the classification. Read pairs with conflicting reads (tag CF) or pairs containing both tags G1 and G2 are considered conflicting and are not reported by default. Reporting of these reads can be enabled using the option --conflicting.
Singleton alignments in the allele-tagged paired-end file (which is the default for e.g. TopHat) are also sorted into the above four files. Specifying --singletons will write these alignments to special singleton files instead (ending in *_st.bam).
Assumes data processed with HiCUP as input, i.e. the input BAM files are by definition paired-end and Reads 1 and 2 follow each other. Hi-C sorting discriminates several more possible read combinations:
G2-G2
G1-UA
G2-UA
G1-G2
UA-UA
Again, read pairs containing a conflicting reads (tag XX:Z:CF) are not printed out by default, but this may be enabled using the option --conflicting. For an example report please see below.
Alignment files produced by the Spliced Transcripts Alignment to a Reference (STAR) aligner also work well with SNPsplit, however a few steps need to be adhered to make this work:
1) Since SNPsplit only recognises the CIGAR operations M, I, D and N (see above), alignments need to be run in end-to-end mode and not using local alignments (which may result in soft-clipping). This can be accomplished using the option: --alignEndsType EndToEnd
Edit 29 09 2020: As of version 0.4.0, SNPsplit also supports soft-clipped reads. (CIGAR operation: S).
2) SNPsplit requires the MD:Z: field of the BAM alignment to work out mismatches involving masked N positions. Since STAR doesn’t report the MD:Z: field by default it needs to be instructed to do so, e.g.: --outSAMattributes NH HI NM MD
3) To save some time and avoid having to sort the reads by name, STAR can be told to leave R1 and R2 following each other in the BAM file using the option: --outSAMtype BAM Unsorted
DNA or RNA alignment files produced by HISAT2 also work well with SNPsplit if you make sure that HISAT2 doesn’t perform soft-clipping. At the time of writing the current version of HISAT2 (2.0.3-beta) does perform soft-clipping (CIGAR operation: S) even though this is not well documented (or in fact the documentation on Github suggests that the default mode is end-to-end which should not perform any soft-clipping whatsoever). Until the end-to-end mode works as expected users will have to set the penalty for soft-clipping so high that it is effectively not performed (--sp is the option governing the soft-clipping penalty). We suggest adding the following option to the HISAT2 command: --sp 1000,1000
Edit: HISAT2 does now also have an option --no-softclip which should have the same effect.
Edit 29 09 2020: As of version 0.4.0, SNPsplit also supports soft-clipped reads. (CIGAR operation: S).
The other very popular Burrows-Wheeler Aligner (BWA) is unfortunately not compatible with SNPsplit alignment sorting as was disscussed in more detail here. Briefly, the reason for this is that BWA randomly replaces the ambiguity nucleobase N in the reference by either C, A, T or G, thereby rendering an N-masked allele-sorting process impossible.
This mode assumes input data has been processed with the bisulfite mapping tool Bismark. SNPsplit will run a quick check at the start of a run to see if the file provided appears to be a Bismark file, and set the flags --bisulfite and/or --paired automatically. In addition it will perform a quick check to see if a paired-end file appears to have been positionally sorted, and if not will set the --no_sort flag (this data is extracted from the @PG header line). Paired-end (--paired) or bisulfite (--bisulfite) mode can also be set manually. Paired-end mode requires Read 1 and Read 2 of a pair to follow each other in consecutive lines. Please be aware that files should not be name-sorted using Sambamba.
In contrast to the standard mode of using all known SNP positions, SNPs involving C to T transitions may not be used for allele-specific sorting since they might reflect either a SNP or a methylation state. This includes all of the following Reference/SNP combinations:
C/T or T/C for forward strand alignments, and G/A or A/G for reverse strand alignments.
The number of SNP positions that have been skipped because of this bisulfite ambiguity is reported in the report file. These positions may however be used to assign opposing strand alignments since they do not involve C to T transitions directly. For that reason, the bisulfite call processing also extracts the bisulfite strand information from the alignments in addition to the base call at the position involved. For any SNPs involving C positions that are not C to T SNPs both methylation states, i.e. C and T, are allowed to match the C position.
For SNPs which are masked by Ns in the genome no methylation call will be performed, i.e. they will receive a . (dot) in the methylation call string. This means that SNP positions are used for allele-sorting only but do not participate in calling methylation. While this may reduce the number of total methylation calls somewhat it effectively eliminates the problem of using positions with potentially incorrect methylation status.
SNPsplit_genome_preparation is designed to read in a variant call file from the Mouse Genomes Project (e.g. this latest file) and generate new genome versions where the strain SNPs are either incorporated into the new genome (full sequence) or masked by the ambiguity nucleobase N (N-masking).
SNPsplit_genome_preparation may be run in two different modes:
1) The VCF file is read and filtered for high-confidence SNPs in the strain specified with strain
2) The reference genome (given with --reference_genome <genome>) is read into memory, and the filtered high-confidence SNP positions are incorporated either as N-masking (default), or full sequence (option --full_sequence)
1) The VCF file is read and filtered for high-confidence SNPs in the strain specified with --strain <name>
2) The reference genome (given with --reference_genome <genome>) is read into memory, and the filtered high-confidence SNP positions are incorporated as full sequence and optionally as N-masking
3) The VCF file is read one more time and filtered for high-confidence SNPs in strain 2 specified with --strain2 <name>
4) The filtered high-confidence SNP positions of strain 2 are incorporated as full sequence, and optionally as N-masking
5) The SNP information of strain and strain 2 relative to the reference genome build are compared, and a new Ref/SNP annotation is constructed whereby the new Ref/SNP information will be Strain/Strain2 (and no longer the standard reference genome strain Black6 (C57BL/6J))
6) The full genome sequence given with --strain <name> is read into memory, and the high-confidence SNP positions between Strain and Strain2 are incorporated as full sequence, and optionally as N-masking
The resulting .fa files are ready to be indexed with your favourite aligner. Proved and tested aligners include Bowtie2, Tophat, STAR, HISAT2, HiCUP and Bismark. Please note that STAR and HISAT2 may require you to disable soft-clipping, please see above for more details.
Both the SNP filtering and the genome preparation write out report files for record keeping.
This section describes in more detail the process of how high confidence SNPs are extracted from the sample VCF file mgp.v5.merged.snps_all.dbSNP142.vcf.gz. We are first going to paste the start of the VCF file and explain in more detail later the individual steps taken by the SNPsplit genome preparation for the strain CAST_EiJ as an example strain; the PDF version of the user guide had relevant information in the header lines or the variant data itself are marked in dark red, but unfotunately this is not supported in this markdown version, so please look carefully... . This should help you adapt the process to other genomes/VCF files should you wish to do so.
VCF file: mgp.v5.merged.snps_all.dbSNP142.vcf.gz
##fileformat=VCFv4.2
##FILTER=<ID=PASS,Description="All filters passed">
##samtoolsVersion=1.1+htslib-1.1
##bcftools_callVersion=1.1+htslib-1.1
##reference=ftp://ftp-mouse.sanger.ac.uk/ref/GRCm38_68.fa
##source_20141009.1=vcf-annotate(r953) -f +/D=50/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: C3H_HeH)
##source_20141009.1=vcf-annotate(r953) -f +/D=100/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: 129S5SvEvBrd,ZALENDE_EiJ,LEWES_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=200/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: C57BL_10J)
##source_20141009.1=vcf-annotate(r953) -f +/D=250/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: 129P2_OlaHsd,A_J,CAST_EiJ,LP_J,PWK_PhJ,WSB_EiJ,BUB_BnJ,DBA_1J,I_LnJ,MOLF_EiJ,NZB_B1NJ,SEA_GnJ,RF_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=300/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: AKR_J,BALB_cJ,C3H_HeJ,C57BL_6NJ,CBA_J,DBA_2J,C57BR_cdJ,C58_J,NZW_LacJ,C57L_J,KK_HiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=350/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: 129S1_SvImJ,FVB_NJ,NOD_ShiLtJ,NZO_HlLtJ,SPRET_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=400/d=5/q=20/w=2/a=5/ (chromosomes=1-19,X,Y: BTBR_T+_Itpr3tf_J,ST_bJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=350/d=5/q=20/w=2/a=5/ (chromosome=MT: LEWES_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=650/d=5/q=20/w=2/a=5/ (chromosome=MT: I_LnJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=850/d=5/q=20/w=2/a=5/ (chromosome=MT: BUB_BnJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=1000/d=5/q=20/w=2/a=5/ (chromosome=MT: SEA_GnJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=1200/d=5/q=20/w=2/a=5/ (chromosome=MT: C57BL_10J)
##source_20141009.1=vcf-annotate(r953) -f +/D=1300/d=5/q=20/w=2/a=5/ (chromosome=MT: RF_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=1450/d=5/q=20/w=2/a=5/ (chromosome=MT: KK_HiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=1550/d=5/q=20/w=2/a=5/ (chromosome=MT: NZB_B1NJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=1800/d=5/q=20/w=2/a=5/ (chromosome=MT: ZALENDE_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=2200/d=5/q=20/w=2/a=5/ (chromosome=MT: C3H_HeH)
##source_20141009.1=vcf-annotate(r953) -f +/D=2300/d=5/q=20/w=2/a=5/ (chromosome=MT: ST_bJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=3050/d=5/q=20/w=2/a=5/ (chromosome=MT: C57L_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=3650/d=5/q=20/w=2/a=5/ (chromosome=MT: MOLF_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=4250/d=5/q=20/w=2/a=5/ (chromosome=MT: NZW_LacJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=4650/d=5/q=20/w=2/a=5/ (chromosome=MT: C3H_HeJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=5300/d=5/q=20/w=2/a=5/ (chromosome=MT: C57BR_cdJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=6150/d=5/q=20/w=2/a=5/ (chromosome=MT: DBA_1J)
##source_20141009.1=vcf-annotate(r953) -f +/D=6200/d=5/q=20/w=2/a=5/ (chromosome=MT: 129S5SvEvBrd,C58_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=6650/d=5/q=20/w=2/a=5/ (chromosome=MT: C57BL_6NJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=6700/d=5/q=20/w=2/a=5/ (chromosome=MT: WSB_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=6900/d=5/q=20/w=2/a=5/ (chromosome=MT: CBA_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=7100/d=5/q=20/w=2/a=5/ (chromosome=MT: BALB_cJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=7450/d=5/q=20/w=2/a=5/ (chromosome=MT: NZO_HlLtJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=7650/d=5/q=20/w=2/a=5/ (chromosome=MT: PWK_PhJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=8300/d=5/q=20/w=2/a=5/ (chromosome=MT: SPRET_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=8850/d=5/q=20/w=2/a=5/ (chromosome=MT: CAST_EiJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=9550/d=5/q=20/w=2/a=5/ (chromosome=MT: 129S1_SvImJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=9600/d=5/q=20/w=2/a=5/ (chromosome=MT: LP_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=9850/d=5/q=20/w=2/a=5/ (chromosome=MT: DBA_2J)
##source_20141009.1=vcf-annotate(r953) -f +/D=10200/d=5/q=20/w=2/a=5/ (chromosome=MT: BTBR_T__Itpr3tf_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=11300/d=5/q=20/w=2/a=5/ (chromosome=MT: AKR_J)
##source_20141009.1=vcf-annotate(r953) -f +/D=11550/d=5/q=20/w=2/a=5/ (chromosome=MT: NOD_ShiLtJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=11700/d=5/q=20/w=2/a=5/ (chromosome=MT: 129P2_OlaHsd)
##source_20141009.1=vcf-annotate(r953) -f +/D=11750/d=5/q=20/w=2/a=5/ (chromosome=MT: FVB_NJ)
##source_20141009.1=vcf-annotate(r953) -f +/D=11800/d=5/q=20/w=2/a=5/ (chromosome=MT: A_J)
##contig=<ID=1,length=195471971>
##contig=<ID=10,length=130694993>
##contig=<ID=11,length=122082543>
##contig=<ID=12,length=120129022>
##contig=<ID=13,length=120421639>
##contig=<ID=14,length=124902244>
##contig=<ID=15,length=104043685>
##contig=<ID=16,length=98207768>
##contig=<ID=17,length=94987271>
##contig=<ID=18,length=90702639>
##contig=<ID=19,length=61431566>
##contig=<ID=2,length=182113224>
##contig=<ID=3,length=160039680>
##contig=<ID=4,length=156508116>
##contig=<ID=5,length=151834684>
##contig=<ID=6,length=149736546>
##contig=<ID=7,length=145441459>
##contig=<ID=8,length=129401213>
##contig=<ID=9,length=124595110>
##contig=<ID=X,length=171031299>
##contig=<ID=Y,length=91744698>
##contig=<ID=MT,length=16299>
##ALT=<ID=X,Description="Represents allele(s) other than observed.">
##QUAL=<ID=QUAL,Number=1,Type=Float,Description="The highest QUAL value for a variant location from any of the samples">
##INFO=<ID=INDEL,Number=0,Type=Flag,Description="Indicates that the variant is an INDEL.">
##INFO=<ID=DP,Number=1,Type=Integer,Description="Raw read depth">
##INFO=<ID=DP4,Number=4,Type=Integer,Description="Total Number of high-quality ref-fwd, ref-reverse, alt-fwd and alt-reverse bases">
##INFO=<ID=CSQ,Number=.,Type=String,Description="Consequence type from Ensembl 75 as predicted by VEP. Format: Allele|Gene|Feature|Feature_type|Consequence|cDNA_position|CDS_position|Protein_position|Amino_acids|Codons|Existing_variation|DISTANCE|STRAND">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=GQ,Number=1,Type=Integer,Description="Phred-scaled Genotype Quality">
##FORMAT=<ID=DP,Number=1,Type=Integer,Description="Number of high-quality bases">
##FORMAT=<ID=MQ0F,Number=1,Type=Float,Description="Fraction of MQ0 reads (smaller is better)">
##FORMAT=<ID=GP,Number=G,Type=Float,Description="Phred-scaled genotype posterior probabilities">
##FORMAT=<ID=PL,Number=G,Type=Integer,Description="List of Phred-scaled genotype likelihoods">
##FORMAT=<ID=AN,Number=1,Type=Integer,Description="Total number of alleles in called genotypes">
##FORMAT=<ID=MQ,Number=1,Type=Integer,Description="Average mapping quality">
##FORMAT=<ID=DV,Number=1,Type=Integer,Description="Number of high-quality non-reference bases">
##FORMAT=<ID=DP4,Number=4,Type=Integer,Description="Number of high-quality ref-fwd, ref-reverse, alt-fwd and alt-reverse bases">
##FORMAT=<ID=SP,Number=1,Type=Integer,Description="Phred-scaled strand bias P-value">
##FORMAT=<ID=SGB,Number=1,Type=Float,Description="Segregation based metric.">
##FORMAT=<ID=PV4,Number=4,Type=Float,Description="P-values for strand bias, baseQ bias, mapQ bias and tail distance bias">
##FORMAT=<ID=FI,Number=1,Type=Integer,Description="Whether a sample was a Pass(1) or fail (0) based on FILTER values">
##FILTER=<ID=StrandBias,Description="Min P-value for strand bias (INFO/PV4) [0.0001]">
##FILTER=<ID=EndDistBias,Description="Min P-value for end distance bias (INFO/PV4) [0.0001]">
##FILTER=<ID=MaxDP,Description="Maximum read depth (INFO/DP or INFO/DP4) []">
##FILTER=<ID=BaseQualBias,Description="Min P-value for baseQ bias (INFO/PV4) [0]">
##FILTER=<ID=MinMQ,Description="Minimum RMS mapping quality for SNPs (INFO/MQ) [20]">
##FILTER=<ID=MinAB,Description="Minimum number of alternate bases (INFO/DP4) [5]">
##FILTER=<ID=Qual,Description="Minimum value of the QUAL field [10]">
##FILTER=<ID=VDB,Description="Minimum Variant Distance Bias (INFO/VDB) [0]">
##FILTER=<ID=GapWin,Description="Window size for filtering adjacent gaps [3]">
##FILTER=<ID=MapQualBias,Description="Min P-value for mapQ bias (INFO/PV4) [0]">
##FILTER=<ID=SnpGap,Description="SNP within INT bp around a gap to be filtered [2]">
##FILTER=<ID=RefN,Description="Reference base is N []">
##FILTER=<ID=MinDP,Description="Minimum read depth (INFO/DP or INFO/DP4) [5]">
##FILTER=<ID=Het,Description="Genotype call is heterozygous (low quality) []">
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT 129P2_OlaHsd 129S1_SvImJ 129S5SvEvBrd AKR_J A_J BALB_cJ BTBR_T+_Itpr3tf_J BUB_BnJ C3H_HeH C3H_HeJ C57BL_10J C57BL_6NJ C57BR_cdJ C57L_J C58_J CAST_EiJ CBA_J DBA_1J DBA_2J FVB_NJ I_LnJ KK_HiJ LEWES_EiJ LP_J MOLF_EiJ NOD_ShiLtJ NZB_B1NJ NZO_HlLtJ NZW_LacJ PWK_PhJ RF_J SEA_GnJ SPRET_EiJ ST_bJ WSB_EiJ ZALENDE_EiJ
1 3000023 . C A 153 MinDP;MinAB;Qual;Het DP=170;DP4=2,0,168,0;CSQ=A||||intergenic_variant|||||||| GT:GQ:DP:MQ0F:GP:PL:AN:MQ:DV:DP4:SP:SGB:PV4:FI 1/1:20:8:0:133,20,0:109,11,0:2:29:7:1,0,7,0:0:-0.6364
26:.:1 1/1:22:6:0.166667:152,22,0:137,18,0:2:36:6:0,0,6,0:0:-0.616816:.:1 1/1:11:4:0:70,11,0:51,4,0:2:24:3:1,0,3,0:0:-0.511536:.:0 1/1:15:4:0.25:80,15,0:68,12,0:2:25:4:0,0,4,0:0:-0.556411:.:0 1/1:11:3:0:72,11,0:63,9,0:2:3
4:3:0,0,3,0:0:-0.511536:.:0 0/1:3:1:0:40,3,3:37,3,0:2:40:1:0,0,1,0:0:-0.379885:.:0 1/1:36:10:0:194,36,0:174,30,0:2:38:10:0,0,10,0:0:-0.670168:.:1 1/1:22:6:0.333333:111,22,0:96,18,0:2:20:6:0,0,6,0:0:-0.616816:.:1 0/1:3:1:0:35,
3,3:32,3,0:2:40:1:0,0,1,0:0:-0.379885:.:0 1/1:19:5:0:135,19,0:121,15,0:2:34:5:0,0,5,0:0:-0.590765:.:1 1/1:15:4:0:106,15,0:94,12,0:2:29:4:0,0,4,0:0:-0.556411:.:0 1/1:11:3:0:83,11,0:74,9,0:2:28:3:0,0,3,0:0:-0.511536:.:0 1/1:33:9:0:199,33,0:180,27,0:2:42:9:0,0,9,0:0:-0.662043:.:1 1/1:6:2:0:56,6,0:50,6,0:2:27:2:0,0,2,0:0:-0.453602:.:0 1/1:15:4:0.5:75,15,0:63,12,0:2:25:4:0,0,4,0:0:-0.556411:.:0 1/1:15:4:0:79,15,0:67,12,0:2:24:4:0,0,4,0:0:-0.556411
:.:0 1/1:22:6:0:133,22,0:118,18,0:2:31:6:0,0,6,0:0:-0.616816:.:1 1/1:15:4:0:108,15,0:96,12,0:2:35:4:0,0,4,0:0:-0.556411:.:0 1/1:15:4:0.25:88,15,0:76,12,0:2:31:4:0,0,4,0:0:-0.556411:.:0 0/0:.:1:0:.,.,.:.,.,.:2:40:1:0,0,1,0:
0:-0.379885:.:0 1/1:11:3:0:79,11,0:70,9,0:2:31:3:0,0,3,0:0:-0.511536:.:0 1/1:6:2:0.5:48,6,0:42,6,0:2:20:2:0,0,2,0:0:-0.453602:.:0 1/1:15:4:0.5:70,15,0:58,12,0:2:22:4:0,0,4,0:0:-0.556411:.:0 1/1:6:2:0:69,6,0:63,6,0:2:40:
2:0,0,2,0:0:-0.453602:.:0 1/1:19:5:0.2:94,19,0:80,15,0:2:24:5:0,0,5,0:0:-0.590765:.:1 1/1:15:4:0:87,15,0:75,12,0:2:30:4:0,0,4,0:0:-0.556411:.:0 1/1:26:7:0.142857:150,26,0:134,21,0:2:31:7:0,0,7,0:0:-0.636426:.:1 1/1:1
9:5:0:132,19,0:118,15,0:2:48:5:0,0,5,0:0:-0.590765:.:1 1/1:30:8:0:171,30,0:153,24,0:2:37:8:0,0,8,0:0:-0.651104:.:1 1/1:11:3:0.333333:72,11,0:63,9,0:2:33:3:0,0,3,0:0:-0.511536:.:0 1/1:19:5:0:137,19,0:123,15,0:2:38:5:0,0,5,0:0:-0.5907
65:.:1 1/1:22:6:0.166667:137,22,0:122,18,0:2:42:6:0,0,6,0:0:-0.616816:.:1 1/1:33:9:0:140,33,0:121,27,0:2:23:9:0,0,9,0:0:-0.662043:.:1 1/1:33:9:0.111111:170,33,0:151,27,0:2:29:9:0,0,9,0:0:-0.662043:.:1 1/1:26:7:0:140,26,0:1
24,21,0:2:31:7:0,0,7,0:0:-0.636426:.:1 1/1:6:2:0:48,6,0:42,6,0:2:21:2:0,0,2,0:0:-0.453602:.:0
1 3000126 rs580370473 G T 184 MinDP;MinMQ;Het;MinAB;Qual DP=417;DP4=93,1,210,113;CSQ=T||||intergenic_variant|||||||| GT:GQ:DP:MQ0F:GP:PL:AN:MQ:DV:DP4:SP:SGB:PV4:FI 1/1:4:8:0:89,4,1:75,0,1:2:36:5:3,0,4,
1:0:-0.590765:.:1 1/1:21:12:0:152,21,0:128,12,0:2:55:9:3,0,7,2:0:-0.662043:.:1 0/0:.:2:0:.,.,.:.,.,.:2:16:0:2,0,0,0:0:.:.:0 0/1:9:11:0:91,0,9:82,0,11:2:48:9:2,0,6,3:0:-0.662043:.:0 1/1:6:5:0:67,6,1:51,1,0:2:44:3:2,0,1,
2:4:-0.511536:.:0 1/1:5:10:0:97,5,1:82,0,0:2:56:8:2,0,7,1:0:-0.651104:.:1 1/1:31:28:0:173,31,0:151,23,0:2:50:21:7,0,6,15:28:-0.692352:.:1 1/1:18:8:0:135,18,0:110,9,0:2:56:6:2,0,3,3:3:-0.616816:.:1 1/1:11:3:0:55,11,0:46,9,0:2:3
7:3:0,0,1,2:0:-0.511536:.:0 1/1:50:14:0:111,50,0:89,42,0:2:41:14:0,0,14,0:0:-0.686358:.:1 1/1:23:7:0:115,23,0:86,12,0:2:38:6:1,0,1,5:0:-0.616816:.:1 1/1:6:8:0:66,6,1:52,2,0:2:49:5:3,0,5,0:0:-0.590765:.:1 1/1:60:17:0.0588235:2
17,60,0:193,51,0:2:52:17:0,0,12,5:0:-0.690438:.:1 1/1:75:33:0:244,75,0:211,62,0:2:52:30:3,0,17,13:6:-0.693097:.:1 1/1:47:13:0:217,47,0:195,39,0:2:49:13:0,0,4,9:0:-0.683931:.:1 0/0:.:5:0:.,.,.:.,.,.:2:47:4:1,0,4,0:0:-0.556411:.:01/1:43:12:0:117,43,0:96,36,0:2:47:12:0,0,12,0:0:-0.680642:.:1 1/1:29:8:0:189,29,0:155,15,0:2:53:7:1,0,3,4:0:-0.636426:.:1 1/1:15:4:0:73,15,0:61,12,0:2:59:4:0,0,4,0:0:-0.556411:.:0 0/1:5:6:0:76,1,5:66,0,7:2:50:4:2,0,4,0:0:-0.5
56411:.:0 1/1:13:9:0:97,13,0:78,7,0:2:45:6:3,0,4,2:3:-0.616816:.:1 1/1:37:15:0.0666667:195,37,0:165,25,0:2:52:13:2,0,7,6:3:-0.683931:.:1 1/1:11:21:0.666667:38,11,0:35,12,0:2:9:6:15,0,3,3:18:-0.616816:.:0 1/1:19:6:0:92
,19,0:67,10,0:2:47:5:1,0,5,0:0:-0.590765:.:1 1/1:10:30:0.4:101,10,0:88,7,0:2:24:14:16,0,0,14:82:-0.686358:.:1 0/1:21:7:0:51,0,21:46,0,21:2:48:5:2,0,4,1:0:-0.590765:.:0 1/1:14:7:0.142857:80,14,0:60,7,0:2:44:5:1,1,4,1:0:-0.
590765:.:1 1/1:5:14:0:101,5,1:86,0,0:2:52:9:5,0,9,0:0:-0.662043:.:1 1/1:33:9:0:196,33,0:177,27,0:2:53:9:0,0,6,3:0:-0.662043:.:1 1/1:26:7:0:91,26,0:75,21,0:2:50:7:0,0,7,0:0:-0.636426:.:1 1/1:44:24:0:138,44,0:119,38,0
:2:45:23:1,0,19,4:0:-0.692717:.:1 1/1:17:8:0:101,17,0:79,9,0:2:53:6:2,0,4,2:0:-0.616816:.:1 0/1:20:11:0.0909091:72,0,20:65,0,21:2:45:7:4,0,7,0:0:-0.636426:.:0 1/1:33:23:0:183,33,0:159,24,0:2:50:18:5,0,6,12:19:-0.691153:.
:1 1/1:11:7:0:66,11,0:49,6,0:2:54:5:2,0,5,0:0:-0.590765:.:1 1/1:19:5:0:62,19,0:48,15,0:2:54:5:0,0,5,0:0:-0.590765:.:1
…
…
To detect the available strains in the VCF file as well as to determine the column number of a desired strain in the file we skim through the header lines until we find a line starting with #CHROM. In terms of the strains in the file the next eight fields are irrelevant:
#CHROM POS ID REF ALT QUAL FILTER INFO FORMAT
The following fields should contain the strains listed in the file, here:
129P2_OlaHsd 129S1_SvImJ 129S5SvEvBrd AKR_J A_J BALB_cJ BTBR_T+_Itpr3tf_J BUB_BnJ C3H_HeH C3H_HeJ C57BL_10J C57BL_6NJ C57BR_cdJ C57L_J C58_J CAST_EiJ CBA_J DBA_1J …
In our example the data for strain CAST_EiJ would be found in column 25, or have an index of 24 (the index = column - 1 because index counting starts at 0 and not at 1).
If #CHROM is not present in the VCF file header, the automated lookup and subsequent steps will fail.
Chromosomes to be processed are detected from the VCF header lines starting with ##contig=<ID=…,… >, e.g. here:
##contig=<ID=1,length=195471971>
##contig=<ID=10,length=130694993>
…
The IDs here are the chromosomes for which variant calls are available, but not necessarily all chromosomes in the genome have to be present here. If there are no variant calls for say “chrY” then the subsequent genome preparation step will simply not introduce any changes but use the reference sequence for “chrY”.
If the ##contig=… fields are missing from the VCF file, subsequent steps are probably going to fail; as a way around this you might get away with defining the chromosome array manually, e.g. like so:
my @chroms = (1..19,’X’,’Y’,’MT’);
For variant calls the SNPsplit genome preparation extracts the following information from each line:
CHROM [Col 1]
POS [Col 2]
REF [Col 4]
ALT [Col 5]
STRAIN [Col 25 (for CAST_EiJ)]
Which in our case is (first three lines only):
CHROM POS REF ALT CAST_EiJ
1 3000023 C A 1/1:15:4:0:79,15,0:67,12,0:2:24:4:0,0,4,0:0:-0.556411:.:0
1 3000126 G T 0/0:.:5:0:.,.,.:.,.,.:2:47:4:1,0,4,0:0:-0.556411:.:0
1 3000185 G T 1/1:43:12:0:276,43,0:255,36,0:2:54:12:0,0,10,2:0:-0.680642:.:1
Now the STRAIN field contains a lot of information which is specified in the FORMAT field and further ‘explained’ in the VCF header. The format is:
GT:GQ:DP:MQ0F:GP:PL:AN:MQ:DV:DP4:SP:SGB:PV4:FI
I am not going into the details about all these FORMAT tags (feel free to browse the header section above), but suffice it to say that the SNPsplit genome preparation only cares about the GT (=GENOTYPE) and FI (=FILTER) entry which are defined as:
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
The GT string can be one of the following:
'.' = no genotype call was made
'0/0' = genotype is the same as the reference genome
'1/1' = homozygous alternative allele; can also be '2/2', '3/3', etc. if more than one alternative allele is present.
'0/1' = heterozygous genotype; can also be '1/2', '0/2', etc.
For the mouse strains we are working with here we only accept homozygous alternative alleles, so '1/1', ‘2/2’ or ‘3/3’. Any other combination is recorded and mentioned in the final report but not included in the SNP file.
##FORMAT=<ID=FI,Number=1,Type=Integer,Description="Whether a sample was a Pass(1) or fail (0) based on FILTER values">
Next we are looking at the FILTER field. If a SNP variant call passed all filters for a given strain it will PASS or get a value of 1. Else it will have a value of 0. In the example above the first position had a genotype of 1/1 but did not pass the Filter for the CAST_EiJ strain (FI value = 0), so this position won’t be included:
1 3000023 C A 1/1:15:4:0:79,15,0:67,12,0:2:24:4:0,0,4,0:0:-0.556411:.:0
The second position did not indicate that there was a SNP (0/0), so the positions wouldn’t be considered for a single-hybrid genome anyway. However, the filter call for that position was uncertain (FI value = 0), so this position would also not be considered for a SNP position in --dual_hybrid mode even if the other strain had a high confidence SNP at this position (this behaviour was introduced in v0.3.2).
1 3000126 G T 0/0:.:5:0:.,.,.:.,.,.:2:47:4:1,0,4,0:0:-0.556411:.:0
The third position finally is a high-confidence homozygous SNP which would be included for both single and dual hybrid genomes:
1 3000185 G T 1/1:43:12:0:276,43,0:255,36,0:2:54:12:0,0,10,2:0:-0.680642:.:1
And finally, variants are required to have a clearly defined genotype, e.g. a made-up position:
12 45630185 A T/C
would not be included as a valid SNP position because the alternative allele is not clearly defined (unless the genotype would be 2/2 here, in which case the reference would be ‘A’ and the alternative allele would be ‘C’).
If you are dealing with any other VCF file than the one used as an example here you might need to adapt one or more of the issues addressed here to make it work.
Input file: 'FVBNJ_Cast.bam'
Writing allele-flagged output file to: 'FVBNJ_Cast.allele_flagged.bam'
Allele-tagging report
=====================
Processed 194564995 read alignments in total
149380724 reads were unassignable (76.78%)
35143075 reads were specific for genome 1 (18.06%)
9860248 reads were specific for genome 2 (5.07%)
118662 reads did not contain one of the expected bases at known SNP positions (0.06%)
180948 contained conflicting allele-specific SNPs (0.09%)
SNP coverage report
===================
N-containing reads: 45276050
non-N: 149262062
total: 194564995
Reads had a deletion of the N-masked position (and were thus dropped): 26883 (0.01%)
Of which had multiple deletions of N-masked positions within the same read: 30
Of valid N containing reads,
N was present in the list of known SNPs: 61087551 (99.99%)
N was not present in the list of SNPs: 4773 (0.01%)
Input file: 'FVBNJ_Cast.allele_flagged.bam'
Writing unassigned reads to: 'FVBNJ_Cast.unassigned.bam'
Writing genome 1-specific reads to: 'FVBNJ_Cast.genome1.bam'
Writing genome 2-specific reads to: 'FVBNJ_Cast.genome2.bam'
Allele-specific paired-end sorting report
=========================================
Read pairs/singletons processed in total: 98215744
thereof were read pairs: 96349251
thereof were singletons: 1866493
Reads were unassignable (not overlapping SNPs): 61174812 (62.29%)
thereof were read pairs: 59662537
thereof were singletons: 1512275
Reads were specific for genome 1: 28657857 (29.18%)
thereof were read pairs: 28446094
thereof were singletons: 211763
Reads were specific for genome 2: 8122687 (8.27%)
thereof were read pairs: 7985424
thereof were singletons: 137263
Reads contained conflicting SNP information: 260388 (0.27%)
thereof were read pairs: 255196
thereof were singletons: 5192
Input file: Black6_129S1.bam
Writing allele-flagged output file to: Black6_129S1.allele_flagged.bam
Allele-tagging report
=====================
Processed 94887256 read alignments in total
59662038 reads were unassignable (62.88%)
19851697 reads were specific for genome 1 (20.92%)
15047281 reads were specific for genome 2 (15.86%)
47261 reads did not contain one of the expected bases at known SNP positions (0.05%)
326240 contained conflicting allele-specific SNPs (0.34%)
SNP coverage report
===================
N-containing reads: 35231977
non-N: 59614777
total: 94887256
Reads had a deletion of the N-masked position (and were thus dropped): 40502 (0.04%)
Of which had multiple deletions of N-masked positions within the same read: 59
Of valid N containing reads,
N was present in the list of known SNPs: 57101748 (99.99%)
N was not present in the list of SNPs: 4211 (0.01%)
Input file: Black6_129S1.allele_flagged.bam'
Writing unassigned reads to: Black6_129S1.UA_UA.bam'
Writing genome 1-specific reads to: Black6_129S1.G1_G1.bam'
Writing genome 2-specific reads to: Black6_129S1.G2_G2.bam'
Writing G1/UA reads to: Black6_129S1.G1_UA.bam'
Writing G2/UA reads to: Black6_129S1.G2_UA.bam'
Writing G1/G2 reads to: Black6_129S1.G1_G2.bam'
Allele-specific paired-end sorting report
=========================================
Read pairs processed in total: 47443628
Read pairs were unassignable (UA/UA): 18862725 (39.76%)
Read pairs were specific for genome 1 (G1/G1): 3533932 (7.45%)
Read pairs were specific for genome 2 (G2/G2): 2592040 (5.46%)
Read pairs were a mix of G1 and UA: 12306421 (25.94%). Of these,
were G1/UA: 6018598
were UA/G1: 6287823
Read pairs were a mix of G2 and UA: 9430675 (19.88%). Of these,
were G2/UA: 4603429
were UA/G2: 4827246
Read pairs were a mix of G1 and G2: 395296 (0.83%). Of these,
were G1/G2: 198330
were G2/G1: 196966
Read pairs contained conflicting SNP information: 322539 (0.68%)
Input file: '129_Cast_bismark_bt2_pe.bam'
Writing allele-flagged output file to: '129_Cast_bismark_bt2_pe.allele_flagged.bam'
Allele-tagging report
=====================
Processed 162441396 read alignments in total
Reads were unaligned and hence skipped: 0 (0.00%)
109109113 reads were unassignable (67.17%)
30267901 reads were specific for genome 1 (18.63%)
22697499 reads were specific for genome 2 (13.97%)
15807753 reads did not contain one of the expected bases at known SNP positions (9.73%)
366883 contained conflicting allele-specific SNPs (0.23%)
SNP coverage report
===================
SNP annotation file: ../all_Cast_SNPs_129S1_reference.mgp.v4.txt.gz
N-containing reads: 68984287
non-N: 93301360
total: 162441396
Reads had a deletion of the N-masked position (and were thus dropped): 155749 (0.10%)
Of which had multiple deletions of N-masked positions within the same read: 65
Of valid N containing reads,
N was present in the list of known SNPs: 119119643 (99.99%)
Positions were skipped since they involved C>T SNPs: 38464451
N was not present in the list of SNPs: 7517 (0.01%)
Input file: 129_Cast_bismark_bt2_pe.allele_flagged.bam'
Writing unassigned reads to: 129_Cast_bismark_bt2_pe.unassigned.bam'
Writing genome 1-specific reads to: 129_Cast_bismark_bt2_pe.genome1.bam'
Writing genome 2-specific reads to: 129_Cast_bismark_bt2_pe.genome2.bam'
Allele-specific paired-end sorting report
=========================================
Read pairs/singletons processed in total: 81220698
thereof were read pairs: 81220698
thereof were singletons: 0
Reads were unassignable (not overlapping SNPs): 40420625 (49.77%)
thereof were read pairs: 40420625
thereof were singletons: 0
Reads were specific for genome 1: 23037433 (28.36%)
thereof were read pairs: 23037433
thereof were singletons: 0
Reads were specific for genome 2: 17303663 (21.30%)
thereof were read pairs: 17303663
thereof were singletons: 0
Reads contained conflicting SNP information: 458977 (0.57%)
thereof were read pairs: 458977
thereof were singletons: 0
USAGE: SNPsplit [options] --snp_file <SNP.file.gz> [input file(s)]
Input file(s)
- Mapping output file in SAM or BAM format. SAM files (ending in
.sam) will first be converted to BAM files.
--snp_file
- Mandatory file specifying SNP positions to be considered, may be a plain text file of gzip compressed. Currently, the SNP file is expected to be in the following format:
SNP-ID Chromosome Position Strand Ref/SNP
33941939 9 68878541 1 T/G
Only the information contained in fields 'Chromosome', 'Position' and 'Ref/SNP base' are being used for analysis. The genome referred to as 'Ref' will be used as genome 1, the genome containing the 'SNP' base as genome 2.
--paired
- Paired-end mode. (Default: OFF).
-o/--outdir <dir>
- Write all output files into this directory. By default the output files will be written into the same folder as the input file(s). If the specified folder does not exist, SNPsplit will attempt to create it first. The path to the output folder can be either relative or absolute.
--singletons
- If the allele-tagged paired-end file also contains singleton alignments (which is the default for e.g. TopHat), these will be written out to extra files (ending in
_st.bam) instead of writing everything to combined paired-end and singleton files. Default: OFF.
--no_sort
- This option skips the sorting step if BAM files are already sorted by read name (e.g. Hi-C files generated by HiCUP). Please note that setting
--no_sortfor unsorted paired-end files will break the tagging process!
--hic
- Assumes Hi-C data processed with HiCUP as input, i.e. the input BAM file is paired-end and Reads 1 and 2 follow each other. Thus, this option also sets the flags
--pairedand--no_sort. Default: OFF.
--bisulfite
- Assumes Bisulfite-Seq data processed with Bismark as input. In paired-end mode (
--paired), Read 1 and Read 2 of a pair are expected to follow each other in consecutive lines. SNPsplit will run a quick check at the start of a run to see if the provided file appears to be a Bismark file, and set the flags--bisulfiteand/or--pairedautomatically. In addition it will perform a quick check to see if a paired-end file appears to have been positionally sorted, and if not will set the flag--no_sort.
--samtools-path
- The path to your Samtools installation, e.g.
/home/user/samtools/. Does not need to be specified explicitly if Samtools is in thePATHenvironment already.
--sam
- The output will be written out in SAM format instead of the default BAM format. SNPsplit will attempt to use the path to Samtools that was specified with
--samtools_path, or, if it hasn't been specified, attempt to find Samtools in thePATHenvironment.
--conflicting/--weird
- Reads or read pairs that were classified as 'Conflicting' (
XX:Z:CF) will be written to an extra file (ending in.conflicting.bam) instead of being simply skipped. Reads may be classified as 'Conflicting' if a single read contains SNP information for both genomes at the same time, or if the SNP position was deleted from the read. Read-pairs are considered 'Conflicting' if either read is was tagged with theXX:Z:CFflag. Default: OFF.
--help
- Displays this help information and exits.
--verbose
(Very!) verbose output (for debugging).
--version
- Displays version information and exits.
USAGE: SNPsplit_genome_preparation [options] --vcf_file <file> --reference_genome /path/to/genome/ --strain <strain_name>
--vcf_file <file>
- Mandatory file specifying SNP information for mouse strains from the Mouse Genomes Project. The file approved and tested is called
mgp.v5.merged.snps_all.dbSNP142.vcf.gz. Please note that future versions of this file or entirely different VCF files might not work out-of-the-box but may require some tweaking. SNP calls are read from the VCF files, and high confidence SNPs are written into a folder in the current working directory calledSNPs_<strain_name>/chr<chromosome>.txt, in the following format:
SNP-ID Chromosome Position Strand Ref/SNP
example: 33941939 9 68878541 1 T/G
--strain <strain_name>
- The strain you would like to use as SNP (ALT) genome. Mandatory. For an overview of strain names just run SNPsplit_genome_preparation selecting
--list_strains.
--list_strains
- Displays a list of strain names present in the VCF file for use with
--strain <strain_name>.
--dual_hybrid
-
Optional. The resulting genome will no longer relate to the original reference specified with
--reference_genome. Instead the new Reference (Ref) is defined by--strain <name>and the new SNP genome is defined by--strain2 <strain_name>.--dual_hybridautomatically sets--full_sequence.This will invoke a multi-step process: 1) Read/filter SNPs for first strain (specified with
--strain < name>) 2) Write full SNP incorporated (and optionally N-masked) genome sequence for first strain 3) Read/filter SNPs for second strain (specified with--strain2 <name>) 4) Write full SNP incorporated (and optionally N-masked) genome sequence for second strain 5) Generate new Ref/Alt SNP annotations for Strain1/Strain2 6) Set first strain as new reference genome and construct full SNP incorporated (and optionally N-masked) genome sequences for Strain1/Strain2
--strain2 <strain_name>
- Optional for constructing dual hybrid genomes (see
--dual_hybridfor more information). For an overview of strain names just run SNPsplit_genome_preparation selecting--list_strains.
--reference_genome
- The path to the reference genome, typically the strain 'Black6' (C57BL/6J), e.g.
--reference_genome /bi/scratch/Genomes/Mouse/GRCm38/. Expects one or more FastA files in this folder (file extension:.faor.fasta).
--skip_filtering
- This option skips reading and filtering the VCF file. This assumes that a folder named
SNPs_<Strain_Name>exists in the working directory, and that text files with SNP information are contained therein in the following format:
SNP-ID Chromosome Position Strand Ref/SNP
example: 33941939 9 68878541 1 T/G
--nmasking
- Write out a genome version for the strain specified where Ref bases are replaced with
N. In the Ref/SNP example T/G the N-masked genome would now carry an N instead of the T. The N-masked genome is written to a folder called<strain_name>_N-masked/. Default: ON.
--full_sequence
- Write out a genome version for the strain specified where Ref bases are replaced with the SNP base. In the Ref/SNP example T/G the full sequence genome would now carry a G instead of the T. The full sequence genome is written out to folder called
<strain_name>_full_sequence/. May be set in addition to--nmasking. Default: OFF.
--no_nmasking
- Disable N-masking if it is not desirable. Will automatically set
--full_sequenceinstead.
--genome_build <name>
- Name of the mouse genome build, e.g. mm10. Will be incorporated into some of the output files. Defaults to
GRCm38.
--help
- Displays this help information and exits.
--version
- Displays version information and exits.