Subcommand: sync

Synopsis: Create a sync file that lists per-sample base counts at each position in the genome.

Usage: grenedalf sync [options]

Documentation for grenedalf v0.4.0

Options

Input SAM/BAM/CRAM

--sam-path

TEXT:PATH(existing)=[] ...
List of sam/bam/cram files or directories to process. For directories, only files with the extension .sam[.gz]|.bam|.cram are processed. To input more than one file or directory, either separate them with spaces, or provide this option multiple times.

--sam-min-map-qual

UINT:UINT in [0 - 90]=0 Needs: --sam-path
Minimum phred-scaled mapping quality score [0-90] for a read in sam/bam/cram files to be considered. Any read that is below the given value of mapping quality will be completely discarded, and its bases not taken into account. Default is 0, meaning no filtering by base quality.

--sam-min-base-qual

UINT:UINT in [0 - 90]=0 Needs: --sam-path
Minimum phred-scaled quality score [0-90] for a base in sam/bam/cram files to be considered. Bases below this are ignored when computing allele frequencies. Default is 0, meaning no filtering by base quality.

--sam-split-by-rg

FLAG Needs: --sam-path
Instead of considering the whole sam/bam/cram file as one large colletion of reads, use the @RG read group tag to split reads. Each read group is then considered a sample. Reads with an invalid (not in the header) read group tag or without a tag are ignored.

--sam-flags-include-all

TEXT Needs: --sam-path
Only use reads with all bits in the given value present in the FLAG field of the read. This is equivalent to the -f / --require-flags setting in samtools view, and uses the same flag names and their corresponding binary values. The value can be specified in hex by beginning with 0x (i.e., /^0x[0-9A-F]+/), in octal by beginning with 0 (i.e., /^0[0-7]+/), as a decimal number not beginning with '0', or as a comma-, plus-, space-, or vertiacal-bar-separated list of flag names as specified by samtools. We are more lenient in parsing flag names than samtools, and allow different capitalization and delimiteres such as dashes and underscores in the flag names as well.

--sam-flags-include-any

TEXT Needs: --sam-path
Only use reads with any bits set in the given value present in the FLAG field of the read. This is equivalent to the --rf / --incl-flags / --include-flags setting in samtools view. See --sam-flags-include-all above for how to specify the value.

--sam-flags-exclude-all

TEXT Needs: --sam-path
Do not use reads with all bits set in the given value present in the FLAG field of the read. This is equivalent to the -G setting in samtools view. See --sam-flags-include-all above for how to specify the value.

--sam-flags-exclude-any

TEXT Needs: --sam-path
Do not use reads with any bits set in the given value present in the FLAG field of the read. This is equivalent to the -F / --excl-flags / --exclude-flags setting in samtools view. See --sam-flags-include-all above for how to specify the value.

Input (m)pileup

--pileup-path

TEXT:PATH(existing)=[] ...
List of (m)pileup files or directories to process. For directories, only files with the extension .(plp|mplp|pileup|mpileup)[.gz] are processed. To input more than one file or directory, either separate them with spaces, or provide this option multiple times.

--pileup-min-base-qual

UINT:UINT in [0 - 90]=0 Needs: --pileup-path
Minimum phred quality score [0-90] for a base in (m)pileup files to be considered. Bases below this are ignored when computing allele frequencies. Default is 0, meaning no filtering by phred quality score.

--pileup-quality-encoding

TEXT:{sanger,illumina-1.3,illumina-1.5,illumina-1.8,solexa}=sanger Needs: --pileup-path
Encoding of the quality scores of the bases in (m)pileup files, when using --pileup-min-base-qual. Default is "sanger", which seems to be the most common these days. Both "sanger" and "illumina-1.8" are identical and use an ASCII offset of 33, while "illumina-1.3" and "illumina-1.5" are identical with an ASCII offset of 64 (we provide different names for completeness). Lastly, "solexa" has an offset of 64, but uses a different equation (not phred score) for the encoding.

Input sync

--sync-path

TEXT:PATH(existing)=[] ...
List of sync (as specified by PoPoolation2) files or directories to process. For directories, only files with the extension .sync[.gz] are processed. To input more than one file or directory, either separate them with spaces, or provide this option multiple times.

Input VCF/BCF

--vcf-path

TEXT:PATH(existing)=[] ...
List of vcf/bcf files or directories to process. For directories, only files with the extension .vcf[.gz]|.bcf are processed. To input more than one file or directory, either separate them with spaces, or provide this option multiple times. This expects that the input file has the per-sample VCF FORMAT field AD (alleleic depth) given, containing the counts of the reference and alternative base. This assumes that the data that was used to create the VCF file was actually a pool of individuals (e.g., from pool sequencing) for each sample (column) of the VCF file. We then interpret the AD field as the allele counts of each pool of individuals. Note that only SNP positions are used; positions that contain indels and other non-SNP variants are skipped.

Input frequency table

--frequency-table-path

TEXT:PATH(existing)=[] ...
List of frequency table files or directories to process. For directories, only files with the extension .(csv|tsv)[.gz] are processed. To input more than one file or directory, either separate them with spaces, or provide this option multiple times.

--frequency-table-separator-char

TEXT:{comma,tab,space,semicolon}=comma Needs: --frequency-table-path
Separator char between fields of the frequency table input.

--frequency-table-depth-factor

FLOAT:POSITIVE=0 Needs: --frequency-table-path
For frequency table input that only contains allele frequencies, without any information on read depth, we need to transform those frequencies into counts for our internal processing. This number is multiplied by the frequency to obtain these pseudo-counts. By default, we use 1000000, to get a reasonable interger approximation of the floating point frequency. This is of course above any typical read depth, but allows for more accurate counts when using for instance haplotype-corrected frequencies such as those from HAF-pipe.

--frequency-table-freq-is-ref

FLAG Needs: --frequency-table-path
For frequency table input that contains allele frequencies, we need to decide whether those frequencies represent the reference or the alternative allele. By default, we assume the latter, i.e., values are interpreted as alternative allele frequencies. Use this flag to instead interpret them as reference allele frequencies.

--frequency-table-chr-column

TEXT Needs: --frequency-table-path
Specify the name of the chromosome column in the header, case sensitive. By default, we look for columns named "chromosome", "chrom", "chr", or "contig", case insensitive.

--frequency-table-pos-column

TEXT Needs: --frequency-table-path
Specify the name of the position column in the header, case sensitive. By default, we look for columns named "position" or "pos", case insensitive.

--frequency-table-ref-base-column

TEXT Needs: --frequency-table-path
Specify the name of the reference base column in the header, case sensitive. By default, we look for columns named "reference", "referencebase", "ref", or "refbase", case insensitive, and ignoring any extra punctuation marks.

--frequency-table-alt-base-column

TEXT Needs: --frequency-table-path
Specify the name of the alternative base column in the header, case sensitive. By default, we look for columns named "alternative", "alternativebase", "alt", or "altbase", case insensitive, and ignoring any extra punctuation marks.

--frequency-table-sample-ref-count-column

TEXT Needs: --frequency-table-path
Specify the exact prefix or suffix of the per-sample reference count columns in the header, case sensitive. By default, we look leniently for column names that combine any of "reference", "referencebase", "ref", or "refbase" with any of "counts", "count", "cnt", or "ct", case insensitive, and ignoring any extra punctuation marks, as a prefix or suffix, with the remainder of the column name used as the sample name. For example, "S1.ref_cnt" indicates the reference count column for sample "S1".

--frequency-table-sample-alt-count-column

TEXT Needs: --frequency-table-path
Specify the exact prefix or suffix of the per-sample alternative count columns in the header, case sensitive. By default, we look leniently for column names that combine any of "alternative", "alternativebase", "alt", or "altbase" with any of "counts", "count", "cnt", or "ct", case insensitive, and ignoring any extra punctuation marks, as a prefix or suffix, with the remainder of the column name used as the sample name. For example, "S1.alt_cnt" indicates the alternative count column for sample "S1".

--frequency-table-sample-freq-column

TEXT Needs: --frequency-table-path
Specify the exact prefix or suffix of the per-sample frequency columns in the header, case sensitive. By default, we look for column names having "frequency", "freq", "maf", "af", or "allelefrequency", case insensitive, and ignoring any extra punctuation marks, as a prefix or suffix, with the remainder of the column name used as the sample name. For example, "S1.freq" indicates the frequency column for sample "S1". Note that when the input data contains frequencies, but no reference or alternative base columns, such as HAF-pipe output tables, we cannot know the bases, and will hence guess. To properly set the reference bases, consider providing the --reference-genome-fasta-file option.

--frequency-table-sample-depth-column

TEXT Needs: --frequency-table-path
Specify the exact prefix or suffix of the per-sample read depth columns in the header, case sensitive. By default, we look for column names having "readdepth", "depth", "coverage", "cov", or "ad", case insensitive, and ignoring any extra punctuation marks, as a prefix or suffix, with the remainder of the column name used as the sample name. For example, "S1.read-depth" indicates the read depth column for sample "S1".

Input Settings

--multi-file-locus-set

TEXT:{union,intersection}=union
When multiple input files are provided, select whether the union of all their loci is used (outer join), or their intersection (inner join). For their union, input files that do not have data at a particular locus are considered as missing at that locus. Note that we allow to use multiple input files even with different file types.

--reference-genome-fasta-file

TEXT:FILE Excludes: --reference-genome-dict-file --reference-genome-fai-file
Provide a reference genome in .fasta[.gz] format. This allows to correctly assign the reference bases in file formats that do not store them, and serves as an integrity check in those that do. It further is used as a sequence dictionary to determine the chromosome order and length, on behalf of a dict or fai file.

--reference-genome-dict-file

TEXT:FILE Excludes: --reference-genome-fasta-file --reference-genome-fai-file
Provide a reference genome sequence dictionary in .dict format. It is used to determine the chromosome order and length, without having to provide the full reference genome.

--reference-genome-fai-file

TEXT:FILE Excludes: --reference-genome-fasta-file --reference-genome-dict-file
Provide a reference genome sequence dictionary in .fai format. It is used to determine the chromosome order and length, without having to provide the full reference genome.

Sample Names, Groups, and Filters

--rename-samples-file

TEXT:FILE
Allows to rename samples, by providing a file that lists the old and new sample names, one per line, separating old and new names by a tab.
By default, we use sample names as provided in the input files. Some file types however do not contain sample names, such as (m)pileup or sync files (unless the non-standard sync header line is provided). For such file types, sample names are automatically assigned by using their input file base name (without path and extension), followed by a dot and numbers 1..n for all samples in that file. For instance, samples in /path/to/sample.sync are named sample.1, sample.2, etc.
Using this option, those names can be renamed as needed. Use verbose output (--verbose) to show a list of all sample names. We then use these names in the output as well as in the --filter-samples-include and --filter-samples-exclude options.

--filter-samples-include

TEXT Excludes: --filter-samples-exclude
Sample names to include (all other samples are excluded); either (1) a comma- or tab-separated list given on the command line (in a typical shell, this list has to be enclosed in quotation marks), or (2) a file with one sample name per line. If no sample filter is provided, all samples in the input file are used. The option is applied after potentially renaming the samples with --rename-samples-file.

--filter-samples-exclude

TEXT Excludes: --filter-samples-include
Sample names to exclude (all other samples are included); either (1) a comma- or tab-separated list given on the command line (in a typical shell, this list has to be enclosed in quotation marks), or (2) a file with one sample name per line. If no sample filter is provided, all samples in the input file are used. The option is applied after potentially renaming the samples with --rename-samples-file.

--sample-group-merge-table-file

TEXT:FILE
When the input contains multiple samples (either within a single input file, or by providing multiple input files), these can be grouped into new samples, by summing up their nucleotide base counts at each position. This has essentially the same effect as having merged the raw fastq files or the mapped sam/bam files of the samples, that is, all reads from those samples are treated as if they were a single sample. For this grouping, the option takes a simple table file (comma- or tab-separated), with the sample names (potentially after the above renaming) in the first column, and their assigned group names in the second column. All samples in the same group are then merged into a grouped sample, and the group names are used as the new sample names for the output. Note that the --pool-sizes option then need to contain the summed up pool sizes for each group.

Region Filters

--filter-region

TEXT=[] ...
Genomic region to filter for, in the format "chr" (for whole chromosomes), "chr:position", "chr:start-end", or "chr:start..end". Positions are 1-based and inclusive (closed intervals). The filter keeps all listed positions, and removes all that are not listed. Multiple region options can be provided, see also --filter-region-set.

--filter-region-list

TEXT:FILE=[] ...
Genomic regions to filter for, as a file with one region per line, either in the format "chr" (for whole chromosomes), "chr:position", "chr:start-end", "chr:start..end", or tab- or space-delimited "chr position" or "chr start end". Positions are 1-based and inclusive (closed intervals). The filter keeps all listed positions, and removes all that are not listed. Multiple region options can be provided, see also --filter-region-set.

--filter-region-bed

TEXT:FILE=[] ...
Genomic regions to filter for, as a BED file. This only uses the chromosome, as well as start and end information per line, and ignores everything else in the file. Note that BED uses 0-based positions, and a half-open [) interval for the end position; simply using columns extracted from other file formats (such as vcf or gff) will not work. The filter keeps all listed positions, and removes all that are not listed.

--filter-region-gff

TEXT:FILE=[] ...
Genomic regions to filter for, as a GFF2/GFF3/GTF file. This only uses the chromosome, as well as start and end information per line, and ignores everything else in the file. The filter keeps all listed positions, and removes all that are not listed.

--filter-region-map-bim

TEXT:FILE=[] ...
Genomic positions to filter for, as a MAP or BIM file as used in PLINK. This only uses the chromosome and coordinate per line, and ignores everything else in the file. The filter keeps all listed positions, and removes all that are not listed.

--filter-region-vcf

TEXT:FILE=[] ...
Genomic positions to filter for, as a VCF/BCF file (such as a known-variants file). This only uses the chromosome and position per line, and ignores everything else in the file. The filter keeps all listed positions, and removes all that are not listed.

--filter-region-fasta

TEXT:FILE=[] ...
Genomic positions to filter for, as a FASTA-like mask file (such as used by vcftools). The file contains a sequence of integer digits [0-9], one for each position on the chromosomes, which specify if the position should be filtered out or not. Any positions with digits above the --filter-region-fasta-min value are removed. Note that this conceptually differs from a mask file, and merely uses the same format.

--filter-region-fasta-min

UINT:INT in [0 - 9]=0 Needs: --filter-region-fasta
When using --filter-region-mask-fasta, set the cutoff threshold for the filtered digits. Only positions with that value or lower will be kept. The default is 0, meaning that all positions with digits greater than 0 will be removed.

--filter-region-fasta-invert

:INT in [0 - 9] Needs: --filter-region-fasta
When using --filter-region-mask-fasta, invert the mask. This option has the same effect as the equivalent in vcftools, but instead of specifying the file, this here is a flag. When it is set, the mask specified above is inverted.

--filter-region-set

TEXT:{union,intersection}=union
It is possible to provide multiple of the above region filter options, even of different types. In that case, decide on how to combine the loci of these filters.

Masking Filters

--filter-mask-bed

TEXT:FILE Excludes: --filter-mask-fasta
Genomic positions to mask (mark as missing), as a BED file. This only uses the chromosome, as well as start and end information per line, and ignores everything else in the file. Note that BED uses 0-based positions, and a half-open [) interval for the end position; simply using columns extracted from other file formats (such as vcf or gff) will not work. The regions listed in the BED file are masked; this is in line with, e.g., smcpp, but differs from the above usage of a BED file for selection regions.

--filter-mask-fasta

TEXT:FILE Excludes: --filter-mask-bed
Genomic positions to mask, as a FASTA-like mask file (such as used by vcftools). The file contains a sequence of integer digits [0-9], one for each position on the chromosomes, which specify if the position should be masked or not. Any positions with digits above the --filter-mask-fasta-min value are tagged as being masked.Note that this conceptually differs from the region fasta above. We here do not remove the masked positions, so that they contribute to, e.g., denominators in the statistics.

--filter-mask-fasta-min

UINT:INT in [0 - 9]=0 Needs: --filter-mask-fasta
When using --filter-mask-fasta, set the cutoff threshold for the masked digits. All positions above that value are masked. The default is 0, meaning that only exactly the positons with value 0 will not be masked.

--filter-mask-fasta-invert

:INT in [0 - 9] Needs: --filter-mask-fasta
When using --filter-mask-fasta, invert the mask. This option has the same effect as the equivalent in vcftools, but instead of specifying the file, this here is a flag. When it is set, the mask specified above is inverted.

Settings

--no-header

FLAG
We provide an extension of the sync format that allows to store sample names in sync files, where a header line is added to the output file of the form: #chr pos ref S1..., where S1... is the list of sample names. Not all other tools that read sync files will be able to parse this, and it hence can be deactivated with this option.

--no-missing-marker

FLAG
We provide an extension of the sync format that allows to mark positions as missing, by using the format .:.:.:.:.:., instead of setting zero counts for these positions. This is particularly useful when storing data from multiple samples, to for instance distinguish true missing data from positions that have been filtered out. Not all other tools that read sync files will be able to parse this, and it hence can be deactivated here, in which case zero counts are written at these positions instead.

--gapless-gsync

FLAG
By default, only the positions for which there is data are printed in the output. However, it might make processing with other tools easier if all files contain all positions, which one might call a gsync file (following the gvcf format). With this option, all missing positions are filled with the missing data indicator, or with zero counts, depending on the --no-missing-marker option. If a referene genome or dictionary is given, this might also include positions beyond where there is input data, up until the length of each chromosome. Hence, the resulting gsync files shall all have the exact same number of lines, which is convenient for simply downstream scripts.

--guess-reference-base

FLAG
By default, when reading from input file formats that do not store the reference base, we do not attempt to guess it. When set however, we use the base with the highest count as the reference base for the output. Alternatively, when a reference genome is provided, we use that to correctly set the reference bases, independently of whether this flag is set.

Output

--out-dir

TEXT=.
Directory to write files to

--file-prefix

TEXT
File prefix for output files. Most grenedalf commands use the command name as the base name for file output. This option amends the base name, to distinguish runs with different data.

--file-suffix

TEXT
File suffix for output files. Most grenedalf commands use the command name as the base name for file output. This option amends the base name, to distinguish runs with different data.

--compress

FLAG
If set, compress the output files using gzip. Output file extensions are automatically extended by .gz.

Global Options

--allow-file-overwriting

FLAG
Allow to overwrite existing output files instead of aborting the command.

--verbose

FLAG
Produce more verbose output.

--threads

UINT
Number of threads to use for calculations. If not set, we guess a reasonable number of threads, by looking at the environmental variables (1) OMP_NUM_THREADS (OpenMP) and (2) SLURM_CPUS_PER_TASK (slurm), as well as (3) the hardware concurrency, taking hyperthreads into account, in the given order of precedence.

--log-file

TEXT
Write all output to a log file, in addition to standard output to the terminal.

Citation

When using this method, please do not forget to cite

Lucas Czech, Jeffrey Spence, Moises Exposito-Alonso. grenedalf: population genetic statistics for the next generation of pool sequencing. arXiv, 2023. doi:10.48550/arXiv.2306.11622

Robert Kofler, Ram Vinay Pandey, Christian Schlötterer. PoPoolation2: identifying differentiation between populations using sequencing of pooled DNA samples (Pool-Seq). Bioinformatics, vol. 27, no. 24, 2011. doi:10.1093/bioinformatics/btr589