You can subscribe to this list here.
2009 |
Jan
|
Feb
(26) |
Mar
(37) |
Apr
(20) |
May
(20) |
Jun
(81) |
Jul
(65) |
Aug
(77) |
Sep
(85) |
Oct
(86) |
Nov
(233) |
Dec
(131) |
---|---|---|---|---|---|---|---|---|---|---|---|---|
2010 |
Jan
(201) |
Feb
(274) |
Mar
(257) |
Apr
(208) |
May
(99) |
Jun
(223) |
Jul
(131) |
Aug
(139) |
Sep
(195) |
Oct
(133) |
Nov
(196) |
Dec
(184) |
2011 |
Jan
(156) |
Feb
(181) |
Mar
(246) |
Apr
(131) |
May
(172) |
Jun
(157) |
Jul
(121) |
Aug
(200) |
Sep
(199) |
Oct
(114) |
Nov
(189) |
Dec
(68) |
2012 |
Jan
(130) |
Feb
(180) |
Mar
(232) |
Apr
(173) |
May
(149) |
Jun
(105) |
Jul
(107) |
Aug
(69) |
Sep
(59) |
Oct
(49) |
Nov
(43) |
Dec
(56) |
2013 |
Jan
(69) |
Feb
(58) |
Mar
(79) |
Apr
(90) |
May
(80) |
Jun
(62) |
Jul
(76) |
Aug
(63) |
Sep
(49) |
Oct
(66) |
Nov
(66) |
Dec
(32) |
2014 |
Jan
(57) |
Feb
(29) |
Mar
(84) |
Apr
(28) |
May
(73) |
Jun
(57) |
Jul
(65) |
Aug
(91) |
Sep
(118) |
Oct
(133) |
Nov
(99) |
Dec
(69) |
2015 |
Jan
(39) |
Feb
(52) |
Mar
(35) |
Apr
(19) |
May
(56) |
Jun
(29) |
Jul
(37) |
Aug
(44) |
Sep
(21) |
Oct
(17) |
Nov
(43) |
Dec
(35) |
2016 |
Jan
(24) |
Feb
(35) |
Mar
(32) |
Apr
(27) |
May
(33) |
Jun
(26) |
Jul
(20) |
Aug
(35) |
Sep
(10) |
Oct
(26) |
Nov
(22) |
Dec
(12) |
2017 |
Jan
(31) |
Feb
(20) |
Mar
(24) |
Apr
(19) |
May
(17) |
Jun
(8) |
Jul
(9) |
Aug
(15) |
Sep
(8) |
Oct
(4) |
Nov
(21) |
Dec
(8) |
2018 |
Jan
(6) |
Feb
(14) |
Mar
(10) |
Apr
(14) |
May
(10) |
Jun
(5) |
Jul
(8) |
Aug
(12) |
Sep
(7) |
Oct
(8) |
Nov
(19) |
Dec
(2) |
2019 |
Jan
(2) |
Feb
(8) |
Mar
(7) |
Apr
|
May
(5) |
Jun
(8) |
Jul
(2) |
Aug
|
Sep
|
Oct
|
Nov
|
Dec
(1) |
2020 |
Jan
|
Feb
(2) |
Mar
(2) |
Apr
(7) |
May
|
Jun
(14) |
Jul
(7) |
Aug
(2) |
Sep
(14) |
Oct
(4) |
Nov
(1) |
Dec
|
2021 |
Jan
|
Feb
|
Mar
(4) |
Apr
(1) |
May
|
Jun
(3) |
Jul
(3) |
Aug
(3) |
Sep
(8) |
Oct
(1) |
Nov
(4) |
Dec
(1) |
2022 |
Jan
(11) |
Feb
(1) |
Mar
|
Apr
(3) |
May
|
Jun
|
Jul
|
Aug
(5) |
Sep
(5) |
Oct
(2) |
Nov
(3) |
Dec
|
2023 |
Jan
|
Feb
(1) |
Mar
|
Apr
|
May
|
Jun
|
Jul
(3) |
Aug
|
Sep
|
Oct
(2) |
Nov
|
Dec
(1) |
2024 |
Jan
(2) |
Feb
|
Mar
|
Apr
(1) |
May
|
Jun
|
Jul
|
Aug
|
Sep
(1) |
Oct
|
Nov
|
Dec
|
From: Robert D. <rm...@sa...> - 2024-09-12 15:59:44
|
Samtools (and HTSlib and BCFtools) version 1.21 is now available from GitHub and SourceForge. https://github.com/samtools/htslib/releases/tag/1.21 https://github.com/samtools/samtools/releases/tag/1.21 https://github.com/samtools/bcftools/releases/tag/1.21 https://sourceforge.net/projects/samtools/ The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.21 ------------------------------------------------------------------------------ The primary user-visible changes in this release are updates to the annot-tsv tool and some speed improvements. Full details of other changes and bugs fixed are below. Notice: this is the last SAMtools / HTSlib release where CRAM 3.0 will be the default CRAM version. From the next we will change to CRAM 3.1 unless the version is explicitly specified, for example using "samtools view -O cram,version=3.0". Updates ------- * Extend annot-tsv with several new command line options. --delim permits use of other delimiters. --headers for selection of other header formats. --no-header-idx to suppress column index numbers in header. Also removed -h as it is now short for --headers. Note --help still works. (PR #1779) * Allow annot-tsv -a to rename annotations. (PR #1709) * Extend annot-tsv --overlap to be able to specify the overlap fraction separately for source and target. (PR #1811) * Added new APIs to facilitate low-level CRAM container manipulations, used by the new "samtools cat" region filtering code. Functions are: cram_container_get_coords() cram_filter_container() cram_index_extents() cram_container_num2offset() cram_container_offset2num() cram_num_containers() cram_num_containers_between() Also improved cram_index_query() to cope with HTS_IDX_NOCOOR regions. (PR #1771) * Bgzip now retains file modification and access times when compressing and decompressing. (PR #1727, fixes #1718. Requested by Gert Hulselmans.) * Use FNV1a for string hashing in khash. The old algorithm was particularly weak with base-64 style strings and lead to a large number of collisions. (PR #1806. Fixes samtools/samtools#2066, reported by Hans-Joachim Ruscheweyh) * Improve the speed of the nibble2base() function on Intel (PR #1667, PR #1764, PR #1786, PR #1802, thanks to Ruben Vorderman) and ARM (PR #1795, thanks to John Marshall). * bgzf_getline() will now warn if it encounters UTF-16 data. (PR #1487, thanks to John Marshall) * Speed up bgzf_read(). While this does not reduce CPU significantly, it does increase the maximum parallelism available permitting 10-15% faster decoding. (PR #1772, PR #1800, Issue #1798) * Speed up faidx by use of better isgraph methods (PR #1797) and whole-line reading (PR #1799, thanks to John Marshall). * Speed up kputll() function, speeding up BAM -> SAM conversion by about 5% and also samtools depth. (PR #1805) * Added more example code, covering fasta/fastq indexing, tabix indexing and use of the thread pool. (PR #1666) Build Changes ------------- * Code warning fixes for pedantic compilers (PR #1777) and avoid some undefined behaviour (PR #1810, PR #1816, PR #1828). * Windows based CI has been migrated from AppVeyor to GitHub Actions. (PR #1796, PR #1803, PR #1808) * Miscellaneous minor build infrastructure and code fixes. (PR #1807, PR #1829, both thanks to John Marshall) * Updated htscodecs submodule to version 1.6.1 (PR #1828) * Fixed an awk script in the Makefile that only worked with gawk. (PR #1831) Bug fixes --------- * Fix small OSS-Fuzz reported issues with CRAM encoding and long CIGARS and/or illegal positions. (PR #1775, PR #1801, PR #1817) * Fix issues with on-the-fly indexing of VCF/BCF (bcftools --write-index) when not using multiple threads. (PR #1837. Fixes samtools/bcftools#2267, reported by Giulio Genovese) * Stricter limits on POS / MPOS / TLEN in sam_parse1(). This fixes a signed overflow reported by OSS-Fuzz and should help prevent other as-yet undetected bugs. (PR #1812) * Check that the underlying file open worked for preload: URLs. Fixes a NULL pointer dereference reported by OSS-Fuzz. (PR #1821) * Fix an infinite loop in hts_itr_query() when given extremely large positions which cause integer overflow. Also adds hts_bin_maxpos() and hts_idx_maxpos() functions. (PR #1774, thanks to John Marshall and reported by Jesus Alberto Munoz Mesa) * Fix an out of bounds read in hts_itr_multi_next() when switching chromosomes. This bug is present in releases 1.11 to 1.20. (PR #1788. Fixes samtools/samtools#2063, reported by acorvelo) * Work around parsing problems with colons in CHROM names. Fixes samtools/bcftools#2139. (PR #1781, John Marshall / James Bonfield) * Correct the CPU detection for Mac OS X 10.7. cpuid is used by htscodecs (see samtools/htscodecs#116), and the corresponding changes in htslib are PR #1785. Reported by Ryan Carsten Schmidt. * Make BAM zero-length intervals work the same as CRAM; permitted and returning overlapping records. (PR #1787. Fixes samtools/samtools#2060, reported by acorvelo) * Replace assert() with abort() in BCF synced reader. This is not an ideal solution, but it gives consistent behaviour when compiling with or without NDEBUG. (PR #1791, thanks to Martin Pollard) * Fixed failure to change the write block size on compressed SAM or VCF files due to an internal type confusion. (PR #1826) * Fixed an out-of-bounds read in cram_codec_iter_next() (PR #1832) ------------------------------------------------------------------------------ samtools - changes v1.21 ------------------------------------------------------------------------------ Notice: * This is the last SAMtools / HTSlib release where CRAM 3.0 will be the default CRAM version. From the next we will change to CRAM 3.1 unless the version is explicitly specified, for example using "samtools view -O cram,version=3.0". New work and changes: * `samtools reset` now removes a set of predefined auxtags, as these tags are no longer valid after the reset operation. This behaviour can be overridden if desired. (PR #2034, fixes #2011. Reported by Felix Lenner) * `samtools reset` now also removes duplicate flags. (PR #2047. Reported by Kevin Lewis) * Region and section/part filtering added to CRAM `samtools cat`. Region filtering permits `samtools cat` to produce new CRAMs that only cover a specified region. (PR #2035) * Added a report of the number of alignments for each primer to `samtools ampliconclip`. (PR #2039, PR #2101, feature request #2033. Thanks to Brad Langhorst) * Make `ampliconclip` primer counts output deterministic. (PR #2081) * `samtools fixmate` no longer removes the PAIRED flag from reads that have no mate. This is done on the understanding that the PAIRED flag is a sequencing technology indicator not a feature of alignment. This is a change to previous `fixmate` behaviour. (PR #2056, fixes #2052. Reported by John Wiedenhoeft) * Added bgzf compressed FASTA output to `samtools faidx`. (PR #2067, fixes #2055. Requested by Filipe G Vieira) * Optimise `samtools depth` histogram incrementing code. (PR #2078) * In `samtools merge` zero pad unique suffix IDs. (PR #2087, fixes #2086. Thanks to Chris Wright) * `samtools idxstats` now accepts the `-X` option, making it easier to specify the location of the index file. (PR #2093, feature request #2071. Requested by Samuel Chen) * Improved documentation for the mpileup `--adjust-MQ` option. (PR #2098. Requested by Georg Langebrake) Bug fixes: * Avoid `tview` buffer overflow for positions with >= 14 digits. (PR #2032. Thanks to John Marshall. Reported on bioconda/bioconda-recipes#47137 by jmunoz94) * Added file name and error message to 'error closing output file' error in `samtools sort`. (PR #2050, fixes #2049. Thanks to Joshua C Randall). * Fixed hard clip trimming issue in `ampliconclip` where right-hand side qualities were being removed from left-hand side trims. (PR #2053, fixes #2048. Reported by Duda5) * Fixed a bug in `samtools merge --template-coordinate` where the wrong heap was being tested. (PR #2062. Thanks to Nils Homer. Reported on ng-core/fastquorum#52 by David Mas-Ponte) * Do not look at chr "*" for unmapped-placed reads with `samtools view --fetch-pairs`. This was causing a significant slowdown when `--fetch-pairs` was being used. (PR #2070, fixes #2059. Reported by acorvelo) * Fixed bug which could cause `samtools view -L` to give incomplete output when the BED file contained nested target locations. (PR #2107, fixes #2104. Reported by geertvandeweyer) * Enable `samtools coverage` to handle alignments that do not have quality score data. This was causing memory access problems. (PR #2083, fixes #2076. Reported by Matthew Colpus) * Fix undefined behaviour in `samtools fastq` with empty QUAL. (PR #2084) * In `plot-bamstats` fixed read-length plot for data with limited variations in length. Lack of data was causing gnuplot problems. (PR #2085, fixes #2068. Reported by mariyeta) * Fixed an accidental fall-through that caused `samtools split -p` to also enable `--no-PG`. (PR #2101) * Fixed an overflow that caused `samtools consensus -m simple` to give incorrect output when the input coverage exceeded several million reads deep. (PR #2099, fixes #2095. Reported by Dylan Lawrence) Non user-visible changes and build improvements: * Work around address sanitizer going missing from the Cirrus CI ubuntu clang compiler by moving the address sanitizer build to gcc. Fix warnings from the new clang compiler. (PR #2043) * Windows based CI has been migrated from AppVeyor to GitHub Actions. (PR #2072, PR #2108) * Turn on more warning options in Cirrus-CI builds, ensure everything builds with `-Werror`, and add undefined behaviour checks to the address sanitizer test. (PR #2101, PR #2103, PR #2109) * Tidy up Makefile dependencies and untracked test files. (PR #2106. Thanks to John Marshall) ------------------------------------------------------------------------------ bcftools - changes v1.21 ------------------------------------------------------------------------------ Changes affecting the whole of bcftools, or multiple commands: * Support multiple semicolon-separated strings when filtering by ID using -i/-e (#2190). For example, `-i 'ID="rs123"'` now correctly matches `rs123;rs456` * The filtering expression ILEN can be positive (insertion), negative (deletion), zero (balanced substitutions), or set to missing value (symbolic alleles). * bcftools query * bcftools +split-vep - The columns indices printed by default with `-H` (e.g., "#[1]CHROM") can be now suppressed by giving the option twice `-HH` (#2152) Changes affecting specific commands: * bcftools annotate - Support dynamic variables read from a tab-delimited annotation file (#2151). For example, in the two cases below the field 'STR' from the -a file is required to match the INFO/TAG in VCF. In the first example the alleles REF,ALT must match, in the second example they are ignored. The option -k is required to output also records that were not annotated: bcftools annotate -a ann.tsv.gz \ -c CHROM,POS,REF,ALT,SCORE,~STR -i'TAG={STR}' -k in.vcf bcftools annotate -a ann.tsv.gz \ -c CHROM,POS,-,-,SCORE,~STR -i'TAG={STR}' -k in.vcf - When adding Type=String annotations from a tab-delimited file, encode characters with special meaning using percent encoding (';', '=' in INFO and ':' in FORMAT) (#2202) * bcftools consensus - Allow to apply a reference allele which overlaps a previous deletion, there is no need to complain about overlapping alleles in such case - Fix a bug which required `-s -` to be present even when there were no samples in the VCF (#2260) * bcftools csq - Fix a rare bug where indel combined with a substitution ending at exon boundary is incorrectly predicted to have 'inframe' rather than 'frameshift' consequence (#2212) * bcftools gtcheck - Fix a segfault with --no-HWE-prob. The bug was introduced with the output format change in 1.19 which replaced the DC section with DCv2 (#2180) - The number of matching genotypes in the DCv2 output was not calculated correctly with non-zero `-E, --error-probability`. Consequently, also the average HWE score was incorrect. The main output, the discordance score, was not affected by the bug * bcftools +mendelian2 - Include the number of good cases where at least one of the trio genotypes has an alternate allele (#2204) - Fix the error message which would report the wrong sample when non-existent sample is given. Note that bug only affected the error message, the program otherwise assigns the family members correctly (#2242) * bcftools merge - Fix a severe bug in merging of FORMAT fields with Number=R and Number=A values. For example, rows with high-coverage FORMAT/AD values (bigger or equal to 128) could have been assigned to incorrect samples. The bug was introduced in version 1.19. For details see #2244. * bcftools mpileup - Return non-zero error code when the input BAM/CRAM file is truncated (#2177) - Add FORMAT/AD annotation by default, disable with `-a -AD` * bcftools norm - Support realignment of symbolic <DUP.*> alleles, similarly to <DEL.*> added previously (#1919,#2145) - Fix in reporting reference allele genotypes with `--multi-overlaps .` (#2160) - Support of duplicate removal of symbolic alleles of the same type but different SVLEN (#2182) - New `-S, --sort` switch to optionally sort output records by allele (#1484) - Add the `-i/-e` filtering options to select records for normalization. Note duplicate removal ignores this option. - Fix a bug where `--atomize` would not fill GT alleles for atomized SNVs followed by an indel (#2239) * bcftools +remove-overlaps - Revamp the program to allow greater flexibility, with the following new options: -M, --mark-tag TAG Mark -m sites with INFO/TAG -m, --mark EXPR Mark (if also -M is present) or remove sites [overlap] dup .. all overlapping sites overlap .. overlapping sites min(QUAL) .. mark sites with lowest QUAL until overlaps are resolved --missing EXPR Value to use for missing tags with -m 'min(QUAL)' 0 .. the default DP .. heuristics, scale maximum QUAL value proportionally to INFO/DP --reverse Apply the reverse logic, for example preserve duplicates instead of removing -O, --output-type t t: plain list of sites (chr,pos), tz: compressed list * bcftools +tag2tag - The conversions --LXX-to-XX, --XX-to-LXX were working but specific cases such as --LAD-to-AD were not. - Print more informative error message when source tag type violiates VCF specification * bcftools +trio-dnm2 - Better handling of the --strictly-novel functionality, especically with respect to chrX inheritance ---------------------------------------------------------------------- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Robert D. <rm...@sa...> - 2024-04-15 17:19:09
|
Samtools (and HTSlib and BCFtools) version 1.20 is now available from GitHub and SourceForge. https://github.com/samtools/htslib/releases/tag/1.20 https://github.com/samtools/samtools/releases/tag/1.20 https://github.com/samtools/bcftools/releases/tag/1.20 https://sourceforge.net/projects/samtools/ The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.20 ------------------------------------------------------------------------------ Updates ------- * When working on named files, bgzip now sets the modified and access times of the output files it makes to match those of the corresponding input. (PR #1727, feature request #1718. Requested by Gert Hulselmans) * It's now possible to use a -o option to specify the output file name in bgzip. (PR #1747, feature request #1726. Requested by Gert Hulselmans) * Improved error faidx error messages. (PR #1743, thanks to Nick Moore) * Faster reading of SAM array (type "B") tags. These often turn up in ONT and PacBio data. (PR #1741) * Improved validity checking of base modification tags. (PR #1749) * mpileup overlap removal now works where one read has a deletion. (PR #1751, fixes samtools/samtools#1992. Reported by Long Tian) * The S3 plugin can now find buckets via S3 access point aliases. (PR #1756, thanks to Matt Pawelczyk; fixes samtools/samtools#1984. Reported by Albert Li) * Added a --threads option (and -@ short option) to tabix. (PR #1755, feature request #1735. Requested by Dan Bolser) * tabix can now index Graph Alignment Format (GAF) files. (See https://github.com/lh3/gfatools/blob/master/doc/rGFA.md) (PR #1763, thanks to Adam Novak) Bug fixes --------- * Security fix: Prevent possible heap overflow in cram_encode_aux() on bad RG:Z tags. (PR #1737) * Security fix: Prevent attempts to call a NULL pointer if certain URL schemes are used in CRAM @SQ UR: tags. (PR #1757) * Security fix: Fixed a bug where following certain AWS S3 redirects could downgrade the connection from TLS (i.e. https://) to unencrypted http://. This could happen when using path-based URLs and AWS_DEFAULT_REGION was set to a region other that the one where the data was stored. (PR #1762, fixes #1760. Reported by andaca) * Fixed arithmetic overflow when loading very long references for CRAM. (PR #1738, fixes #1738. Reported by Shane McCarthy) * Fixed faidx and CRAM reference look-ups on compressed fasta where the .fai index file was present, but the .gzi index of compressed offsets was not. (PR #1745, fixes #1744. Reported by Theodore Li) * Fixed BCF indexing on-the-fly bug which produced invalid indexes when using multiple compression threads. (PR #1742, fixes #1740. Reported by graphenn) * Ensure that pileup destructors are called by bam_plp_destroy(), to prevent memory leaks. (PR #1749, PR #1754) * Ensure on-the-fly index timestamps are always older than the data file. Previously the files could be closed out of order, leading to warnings being printed when using the index. (PR #1753, fixes #1732. Reported by Gert Hulselmans) * To prevent data corruption when reading (strictly invalid) VCF files with duplicated FORMAT tags, all but the first copy of the data associated with the tag are now dropped with a warning. (PR #1752, PR #1761, fixes #1733. Reported by anthakki) * Fixed a bug introduced in release 1.19 (PR #1689) which broke variant record data if it tried to remove an over-long tag. (PR #1752, PR #1761) * Changed error to warning when complaining about use of the CG tag in SAM or CRAM files. (PR #1758, fixes samtools/samtools#2002) ------------------------------------------------------------------------------ samtools - changes v1.20 ------------------------------------------------------------------------------ * Added a `--max-depth` option to `bedcov`, for more control over the depth limit used when calculating the pileup. Previously this limit was set at 64000; now it is set to over 2 billion, so effectively all bases will be counted. (PR #1970, fixes #1950. Reported by ellisjj) * Added `mpileup --output-extra RLEN` to display the unclipped read length. (PR #1971, feature request #1959. Requested by Feng Tian) * Improved checking of symbolic flag names (e.g. UNMAP) passed to samtools. (PR #1981, fixes #1977. Reported by Ilya Shlyakhter) * The `samtools consensus --min-depth` option now works for the Bayesian mode as well as the simple one. (PR #1989, feature request #1982. Requested by Gautier Richard) * It's now possible to use the `samtools fastq` `-d tag:val` option multiple times, allowing matches on more than one tag/value. It also gets a `-D` option which allows the values to be listed in a file. (PR #1993, feature request #1958. Requested by Tristan Lefebure) * Added `samtools fixmate` `-M` option to sanity check base modification (`ML`, `MM`, `MN`) tags, and where necessary adjust modification data on hard-clipped records. (PR #1990) * Made `mpileup` run faster. (PR #1995) * `samtools import` now adds a `@PG` header to the files it makes. As with other sub-commands, this can be disabled by using `--no-PG`. (PR #2008. Requested by Steven Leonard) * The `samtools split` `-d` option to split by tag value now works on tags with integer values. (PR #2005, feature request #1956. Requested by Alex Leonard) * Adjusted `samtools sort -n` (by name) so that primary reads are always sorted before secondary / supplementary. (PR #2012, feature request #2010. Requested by Stijn van Dongen) * Added `samtools bedcov` `-H` option to print column headers in the output. (PR #2025. Thanks to Dr. K. D. Murray) Documentation: * Added a note that BAQ is applied before filtering and overlap removal during mpileup processing. (PR #1988, fixes #1985. Reported by Joseph Galasso) * Added 3.1 to the list of supported CRAM versions in the samtools manual page. (PR #2009. Thanks to Andrew Thrasher) * Made assorted improvements to ampliconclip, flagstat and markdup manual pages. (PR #2014) Bug Fixes: * Security fix: Fixed double free that could occur if bed file indexing failed due to running out of memory. This bug first appeared in version 1.19.1. (PR #2026) * Corrected error message printed when faidx fails to load the fai index. (PR #1987. Thanks to Nick Moore) * Fixed bug introduced in release 1.4 that caused incorrect reference bases to be printed by `samtools mpileup -a -f ref.fa` in the zero-depth regions at the end of each reference. (PR #2019, fixes #2018. Reported by Joe Georgeson) * Fixed a samtools view usage crash on MinGW when given invalid options. (PR #2030, fixes #2029. Reported by Divon Lan) Non user-visible changes and build improvements: * Added tests to ensure that CRAM compression is working properly. (PR #1969, part of fix for #1968. Reported by Clockris) ------------------------------------------------------------------------------ bcftools - changes v1.20 ------------------------------------------------------------------------------ Changes affecting the whole of bcftools, or multiple commands: * Add short option -W for --write-index. The option now accepts an optional parameter which allows to choose between TBI and CSI index format. Changes affecting specific commands: * bcftools consensus - Add new --regions-overlap option which allows to take into account overlapping deletions that start out of the fasta file target region. * bcftools isec - Add new option `-l, --file-list` to read the list of file names from a file * bcftools merge - Add new option `--force-single` to support single-file edge case (#2100) * bcftools mpileup - Add new option --indels-cns for an alternative indel calling model, which should increase the speed on long read data (thanks to using edlib) and the precision (thanks to a number of heuristics). * bcftools norm - Change the order of atomization and multiallelic splitting (when both -a,-m are given) from "atomize first, then split" to "split first, then atomize". This usually results in a simpler VCF representation. The previous behaviour can be achieved by explicitly streaming the output of the --atomize command into the --multiallelics splitting command. - Fix Type=String multiallelic splitting for Number=A,R,G tags with incorrect number of values. - Merging into multiallelic sites with `bcftools norm -m +indels` did not work. This is now fixed and the merging is now more strict about variant types, for example complex events, such as AC>TGA, are not considered as indels anymore (#2084) * bcftools reheader - Allow reading the input file from a stream with --fai (#2088) * bcftools +setGT - Support for custom genotypes based on the allele with higher depth, such as `--new-gt c:0/X` custom genotypes (#2065) * bcftools +split-vep - When only one of the tags is present, automatically choose INFO/BCSQ (the default tag name produced by `bcftools csq`) or INFO/CSQ (produced by VEP). When both tags are present, use the default INFO/CSQ. - Transcript selection by MANE, PICK, and user-defined transcripts, for example --select CANONICAL=YES --select MANE_SELECT!="" --select PolyPhen~probably_damaging - Select all matching transcripts via --select, not just one - Change automatic type parsing of VEP fields DNA_position, CDS_position, and Protein_position from Integer to String, as it can be of the form "8586-8599/9231". The type Integer can be still enforced with `-c cDNA_position:int,CDS_position:int,Protein_position:int`. - Recognize `-c field:str`, not just `-c field:string`, as advertised in the usage page - Fix a bug which made filtering expression containing missing values crash (#2098) * bcftools stats - When GT is missing but AD is present, the program determines the alternate allele from AD. However, if the AD tag has incorrect number of values, the program would exit with an error printing "Requested allele outside valid range". This is now fixed by taking into account the actual number of ALT alleles. * bcftools +tag2tag - Support for conversion from tags using localized alleles (e.g. LPL, LAD) to the family of standard tags (PL, AD) * bcftools +trio-dnm2 - Extend --strictly-novel to exclude cases where the non-Mendelian allele is the reference allele. The change is motivated by the observation that this class of variants is enriched for errors (especially for indels), and better corresponds with the option name. -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Robert D. <rm...@sa...> - 2024-01-24 12:43:15
|
Samtools release 1.19.2 is now available from GitHub and SourceForge. It fixes an error in 1.19.1 that broke filtering on unordered BED files. https://github.com/samtools/samtools/releases/tag/1.19.2 https://sourceforge.net/projects/samtools/ ------------------------------------------------------------------------------ samtools - changes v1.19.2 ------------------------------------------------------------------------------ Bug Fixes: * Fixed a regression in 1.19.1 that broke BED filtering for inputs where the region start positions for the same reference were not sorted in ascending order. (PR #1975, fixes #1974. Reported by Anže Starič) -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Robert D. <rm...@sa...> - 2024-01-22 12:10:43
|
Samtools and HTSlib release 1.19.1 are now available from GitHub and SourceForge. This release fixes a regression in version 1.19 that caused written CRAM files to be much bigger than necessary. It also fixes a number of other bugs listed below, including one that could cause crashes or incorrect results when filtering by regions in a BED file. https://github.com/samtools/htslib/releases/tag/1.19.1 https://github.com/samtools/samtools/releases/tag/1.19.1 https://sourceforge.net/projects/samtools/ ------------------------------------------------------------------------------ htslib - changes v1.19.1 ------------------------------------------------------------------------------ * Fixed a regression in release 1.19 that caused all aux records to be stored uncompressed in CRAM files. The resulting files were correctly formatted, but bigger than they needed to be. (PR#1729, fixes samtools#1968. Reported by Clockris) * Fixed possible out-of-bounds reads due to an incorrect check on B tag lengths in cram_encode_aux(). (PR#1725) * Fixed an incorrect check on tag length which could fail to catch a two byte out-of-bounds read in bam_get_aux(). (PR#1728) * Made errors reported by hts_open_format() less confusing when it can't open the reference file. (PR#1724, fixes #1723. Reported by Alex Leonard) * Made hts_close() fail more gracefully if it's passed a NULL pointer (PR#1724) ------------------------------------------------------------------------------ samtools - changes v1.19.1 ------------------------------------------------------------------------------ Bug Fixes: * Fixed a possible array bounds violation when looking up regions in a BED file (e.g. using `samtools view -L`). This could lead to crashes or the return of incomplete results if the BED file contained a large number of entries all referencing low positions on a chromosome. (PR #1962, fixes #1961. Reported by geertvandeweyer) * Fixed a crash in samtools stats that occurred when trying to clean up after it was unable to open a CRAM reference file. (PR #1957, fixes crash reported in samtools/htslib#1723. Reported by Alex Leonard) Documentation: * Fixed inverted logic in the `samtools consensus --show-del` manual page description. (PR #1955, fixes #1951. Reported by Mikhail Schelkunov) * Added a description of the MPC section to the `samtools stats` manual page. (PR #1963, fixes #1954. Reported by litun-fkby) -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Robert D. <rm...@sa...> - 2023-12-12 16:23:04
|
Samtools (and HTSlib and BCFtools) version 1.19 is now available from GitHub and SourceForge. https://github.com/samtools/htslib/releases/tag/1.19 https://github.com/samtools/samtools/releases/tag/1.19 https://github.com/samtools/bcftools/releases/tag/1.19 https://sourceforge.net/projects/samtools/ The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.19 ------------------------------------------------------------------------------ Updates ------- * A temporary work-around has been put in the VCF parser so that it is less likely to fail on rows with a large number of ALT alleles, where Number=G tags like PL can expand beyond the 2Gb limit enforced by HTSlib. For now, where this happens the offending tag will be dropped so the data can be processed, albeit without the likelihood data. In future work, the library will instead convert such tags into their local alternatives (see https://github.com/samtools/hts-specs/pull/434). * New program. Adds annot-tsv which annotates regions in a destination file with texts from overlapping regions in a source file. (PR#1619) * Change bam_parse_cigar() so that it can modify existing BAM records. This makes more useful as public API. Previously it could only handle partially formed BAM records. (PR#1651, fixes #1650. Reported by Oleksii Nikolaienko) * Add "uncompressed" to hts_format_description() where appropriate. This adds an "uncompressed" description to uncompressed files that would normally be compressed, such as BAM and BCF. (PR#1656, in relation to samtools#1884. Thanks to John Marshall) * Speed up to the VCF parser and writer. (PR#1644 and PR#1663) * Add an hclen (hard clip length) SAM filter function. (PR#1660, with reference to samtools#813) * Avoid really closing stdin/stdout in hclose()/hts_close()/et al. See discussion in PR for details. (PR#1665. Thanks to John Marshall) * Add support to handle multiple files in bgzip. (PR#1658, fixes #1642. Requested by bw2) * Enable auto-vectorisation in CRAM 3.1 codecs. Speeds decoding on some sequencing platform data. (PR#1669) * Speed up removal of lines in large headers. (PR#1662, fixes #1460. Reported by Anže Starič) * Apply seqtk PR to improve kseq.h parsing performance. Port of Fabian Klötzl's (kloetzl) lh3/seqtk#123 and attractivechaos/klib#173 to HTSlib. (PR#1674. Thanks to John Marshall) Build changes ------------- * Updated htscodecs submodule to 1.6.0. (PR#1685, PR#1717, PR#1719) * Apply the packed attribute to uint*_u types for Clang to prevent -fsanitize=alignment failures. (PR#1667. Thanks to Fangrui Song) * Fuzz testing improvements. (PR#1664) * Add C++ casts for external headers in klist.h and kseq.h. (PR#1683. See also PR#1674 and PR#1682) * Add test case compiling the public headers as C++. (PR#1682. Thanks to John Marshall) * Enable optimisation level -O3 for SAM QUAL+33 formatting. (PR#1679) * Make compiler flag detection work with zig cc. (PR#1687) * Fix unused value warnings when built with NDEBUG. (PR#1688) * Remove some disused Makefile variables, fix typos and a warning. Improve bam_parse_basemod() documentation. (PR#1705, Thanks to John Marshall) Bug fixes --------- * Fail bgzf_useek() when offset is above block limits. (PR#1668) * Fix multi-threaded on-the-fly indexing problems. (PR#1672, fixes samtools#1861 and bcftools#1985. Reported by Mark Ebbert and lacek) * Fix hfile_libcurl small seek bug. (PR#1676, fixes samtools#1918. Also may fix #1037, #1625 and samtools#1622. Reported by Alex Reynolds, Mark Walker, Arthur Gilly and skatragadda-nygc. Thanks to John Marshall) * Fix a minor memory leak in malformed CRAM EXTERNAL blocks. [fuzz] (PR#1671) * Fix a cram decode hang from block_resize(). (PR#1680. Reported by Sebastian Deorowicz) * Cram fuzzing improvements. Fixes a number of cram errors. (PR#1701, fixes #1691, #1692, #1693, #1696, #1697, #1698, #1699 and #1700. Thanks to Octavio Galland for finding and reporting all these) * Fix crypt4gh redirection. (PR#1675, fixes grbot/crypt4gh-tutorial#2. Reported by hth4) * Fix PG header linking when records make a loop. (PR#1702, fixes #1694. Reported by Octavio Galland) * Prevent issues with no-stored-sequence records in CRAM files, by ensuring they are accounted for properly in block size calculations, and by limiting the maximum query length in the CIGAR data. Originally seen as an overflow by OSS-Fuzz / UBSAN, it turned out this could lead to excessive time and memory use by HTSlib, and could result in it writing out unreadable CRAM files. (PR#1710) * Fix some illegal shifts and integer overflows found by OSS-Fuzz / UBSAN. (PR#1707, PR#1712, PR#1713) ------------------------------------------------------------------------------ samtools - changes v1.19 ------------------------------------------------------------------------------ New work and changes: * Samtools coverage: add a new --plot-depth option to draw depth (of coverage) rather than the percentage of bases covered. (PR #1910. Thanks to Pierre Lindenbaum) * Samtools merge / sort: add a lexicographical name-sort option via the -N option. The "natural" alpha-numeric sort is still available via -n. (PR #1900, fixes #1500. Reported by Steve Huang) * Samtools view: add -N ^NAME_FILE and -R ^RG_FILE options. The standard -N and -R options only output reads matching a specified list of names or read-groups. With a caret (^) prefix these may be negated to only output read not matching the specified files. (PR #1896, fixes #1895. Suggested by Feng Tian) * Cope with Htslib's change to no longer close stdout on hts_close. Htslib companion PR is samtools/htslib#1665. (PR #1909. Thanks to John Marshall) * Plot-bamstats: add a new plot of the read lengths ("RL") from samtools stats output. (PR #1922, fixes #1824. Thanks to @erboone, suggested by Martin Pollard) * Samtools split: support splitting files based on the contents of auxiliary tags. Also adds a -M option to limit the number of files split can make, to avoid accidental resource over-allocation, and fixes some issues with --write-index. (PR #1222, PR #1933, fixes #1758. Thanks to Valeriu Ohan, suggested by Scott Norton) Bug Fixes: * Samtools stats: empty barcode tags are now treated as having no barcode. (PR #1929, fixes #1926. Reported by Jukka Matilainen) * Samtools cat: add support for non-seekable streams. The file format detection blocked pipes from working before, but now files may be non-seekable such as stdin or a pipe. (PR #1930, fixes #1731. Reported by Julian Hess) * Samtools mpileup -aa (absolutely all positions) now produces an output even when given an empty input file. (PR #1939. Reported by Chang Y) * Samtools markdup: speed up optical duplicate tagging on regions with very deep data. (PR #1952) Documentation: * Samtools mpileup: add more usage examples to the man page. (PR #1913, fixes #1801) * Samtools fastq: explicitly document the order that filters apply. (PR #1907) * Samtools merge: fix example output to use an uppercase RG PL field. (PR #1917. Thanks to John Marshall. Reported by Michael Macias) * Add hclen SAM filter documentation. (PR #1902. See also samtools/htslib#1660) * Samtools consensus: remove the -5 option from documentation. This option was renamed before the consensus subcommand was merged, but accidentally left in the man page. (PR #1924) ------------------------------------------------------------------------------ bcftools - changes v1.19 ------------------------------------------------------------------------------ Changes affecting the whole of bcftools, or multiple commands: * Filtering expressions can be given a file with list of strings to match, this was previously possible only for the ID column. For example ID=@file .. selects lines with ID present in the file INFO/TAG=@file.txt .. selects lines where TAG has a string value listed in the file INFO/TAG!=@file.txt .. TAG must not have a string value listed in the file * Allow to query REF,ALT columns directly, for example -e 'REF="N"' Changes affecting specific commands: * bcftools annotate - Fix `bcftools annotate --mark-sites`, VCF sites overlapping regions in a BED file were not annotated (#1989) - Add flexibility to FILTER column transfers and allow transfers within the same file, across files, and in combination. For examples see http://samtools.github.io/bcftools/howtos/annotate.html#transfer_filter_to_info * bcftools call - Output MIN_DP rather than MinDP in gVCF mode - New `-*, --keep-unseen-allele` option to output the unobserved allele <*>, intended for gVCF. * bcftools head - New `-s, --samples` option to include the #CHROM header line with samples. * bcftools gtcheck - Add output options `-o, --output` and `-O, --output-type` - Add filtering options `-i, --include` and `-e, --exclude` - Rename the short option `-e, --error-probability` from lower case to upper case `-E, --error-probability` - Changes to the output format, replace the DC section with DCv2: - adds a new column for the number of matching genotypes - The --error-probability is newly interpreted as the probability of erroneous allele rather than genotype. In other words, the calculation of the discordance score now considers the probability of genotyping error to be different for HOM and HET genotypes, i.e. P(0/1|dsg=0) > P(1/1|dsg=0). - fixes in HWE score calculation plus output average HWE score rather than absolute HWE score - better description of fields * bcftools merge - Add `-m` modifiers to suppress the output of the unseen allele <*> or <NON_REF> at variant sites (e.g. `-m both,*`) or all sites (e.g. `-m both,**`) * bcftools mpileup - Output MIN_DP rather than MinDP in gVCF mode * bcftools norm - Add the number of joined lines to the summary output, for example Lines total/split/joined/realigned/skipped: 6/0/3/0/0 - Allow combining -m and -a with --old-rec-tag (#2020) - Symbolic <DEL> alleles caused norm to expand REF to the full length of the deletion. This was not intended and problematic for long deletions, the REF allele should list one base only (#2029) * bcftools query - Add new `-N, --disable-automatic-newline` option for pre-1.18 query formatting behavior when newline would not be added when missing - Make the automatic addition of the newline character in a more predictable way and, when missing, always put it at the end of the expression. In version 1.18 it could be added at the end of the expression (for per-site expressions) or inside the square brackets (for per-sample expressions). The new behavior is: - if the formatting expression contains a newline character, do nothing - if there is no newline character and -N, --disable-automatic-newline is given, do nothing - if there is no newline character and -N is not given, insert newline at the end of the expression See #1969 for details - Add new `-F, --print-filtered` option to output a default string for samples that would otherwise be filtered by `-i/-e` expressions. - Include sample name in the output header with `-H` whenever it makes sense (#1992) * bcftools +spit-vep - Fix on the fly filtering involving numeric subfields, e.g. `-i 'MAX_AF<0.001'` (#2039) - Interpret default column type names (--columns-types) as entire strings, rather than substrings to avoid unexpected spurious matches (i.e. internally add ^ and $ to all field names) * bcftools +trio-dnm2 - Do not flag paternal genotyping errors as de novo mutations. Specifically, when father's chrX genotype is 0/1 and mother's 0/0, 0/1 in the child will not be marked as DNM. * bcftools view - Add new `-A, --trim-unseen-allele` option to remove the unseen allele <*> or <NON_REF> at variant sites (`-A`) or all sites (`-AA`) -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Tom B. <tb...@um...> - 2023-10-24 16:09:39
|
Chris - There must be oodles of journal articles over the years since about 2010 in the biology literatrure describing this analysis for non-human species. I have no experience myself with non-human species. A very vanilla workflow using bwa-mem and samtools is briefly described at https://www.htslib.org/workflow/wgs-call.html. The Broad Institute GATK suite of tools has very extensive and detailed documentation. https://gatk.broadinstitute.org/hc/en-us/articles/360035535932-Germline-short-variant-discovery-SNPs-Indels- is one example. But as an introduction, the right journal article will probably be more useful than either of these technical web pages. - thanks - tom blackwell - On Tue, Oct 24, 2023 at 6:10 AM Christopher Blair <cb...@gc...> wrote: > Hi all, > > I was hoping that someone could point me in the right direction for > processing WGS data for phylogenetics and population genetics. I will soon > have 10x whole genome data for about 27 individuals of a lizard species. I > will need to go from raw reads to multiple sequence alignments and > SNPs/genotypes. There is a reference genome available from a closely > related species. Any guidance would be appreciated. > > > Chris > > > > > <{===================~~~~~~~~~~~ > /\ /\ > Christopher Blair, Ph.D. > Executive Officer, Biology PhD Program > > Professor > > Graduate Center and > > New York City College of Technology > The City University of New York > cb...@ci... <cb...@gc...>; cb...@gc... > Website: https://sites.google.com/site/christopherblairphd/home > Office: 4315 (GC); Academic Complex 501D (City Tech) > Lab: Academic Complex 401 (City Tech) > _______________________________________________ > Samtools-help mailing list > Sam...@li... > https://lists.sourceforge.net/lists/listinfo/samtools-help > |
From: Christopher B. <cb...@gc...> - 2023-10-23 19:28:21
|
Hi all, I was hoping that someone could point me in the right direction for processing WGS data for phylogenetics and population genetics. I will soon have 10x whole genome data for about 27 individuals of a lizard species. I will need to go from raw reads to multiple sequence alignments and SNPs/genotypes. There is a reference genome available from a closely related species. Any guidance would be appreciated. Chris <{===================~~~~~~~~~~~ /\ /\ Christopher Blair, Ph.D. Executive Officer, Biology PhD Program Professor Graduate Center and New York City College of Technology The City University of New York cb...@ci...<mailto:cb...@gc...>; cb...@gc...<mailto:cb...@gc...> Website: https://sites.google.com/site/christopherblairphd/home Office: 4315 (GC); Academic Complex 501D (City Tech) Lab: Academic Complex 401 (City Tech) |
From: Robert D. <rm...@sa...> - 2023-07-25 13:14:22
|
Samtools (and HTSlib and BCFtools) version 1.18 is now available from GitHub and SourceForge. https://github.com/samtools/htslib/releases/tag/1.18 https://github.com/samtools/samtools/releases/tag/1.18 https://github.com/samtools/bcftools/releases/tag/1.18 https://sourceforge.net/projects/samtools/ The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.18 ------------------------------------------------------------------------------ Updates ------- * Using CRAM 3.1 no longer gives a warning about the specification being draft. Note CRAM 3.0 is still the default output format. (PR#1583) * Replaced use of sprintf with snprintf, to silence potential warnings from Apple's compilers and those who implement similar checks. (PR#1594, fixes #1586. Reported by Oleksii Nikolaienko) * Fastq output will now generate empty records for reads with no sequence data (i.e. sequence is "*" in SAM format). (PR#1576, fixes samtools/samtools#1576. Reported by Nils Homer) * CRAM decoding speed-ups. (PR#1580) * A new MN aux tag can now be used to verify that MM/ML base modification data has not been broken by hard clipping. (PR#1590, PR#1612. See also PR samtools/hts-specs#714 and issue samtools/hts-specs#646. Reported by Jared Simpson) * The base modification API has been improved to make it easier for callers to tell unchecked bases from unmodified ones. (PR#1636, fixes #1550. Requested by Chris Wright) * A new bam_mods_queryi() API has been added to return additional data about the i-th base modification returned by bam_mods_recorded(). (PR#1636, fixes #1550 and #1635. Requested by Jared Simpson) * Speed up index look-ups for whole-chromosome queries. (PR#1596) * Mpileup now merges adjacent (mis)match CIGAR operations, so CIGARs using the X/= operators give the same results as if the M operator was used. (PR#1607, fixes #1597. Reported by Marcel Martin) * It's now possible to call bcf_sr_set_regions() after adding readers using bcf_sr_add_reader() (previously this returned an error). Doing so will discard any unread data, and reset the readers so they iterate over the new regions. (PR#1624, fixes samtools/bcftools#1918. Reported by Gregg Thomas) * The synced BCF reader can now accept regions with reference names including colons and hyphens, by enclosing them in curly braces. For example, {chr_part:1-1001}:10-20 will return bases 10 to 20 from reference "chr_part:1-1001". (PR#1630, fixes #1620. Reported by Bren) * Add a "samples" directory with code demonstrating usage of HTSlib plus a tutorial document. (PR#1589) Build changes ------------- * Htscodecs has been updated to 1.5.1 (PR#1654) * Htscodecs SIMD code now works with Apple multiarch binaries. (PR#1587, HTSlib fix for samtools/htscodecs#76. Reported by John Marshall) * Improve portability of "expr" usage in version.sh. (PR#1593, fixes #1592. Reported by John Marshall) * Improve portability to *BSD targets by ensuring _XOPEN_SOURCE is defined correctly and that source files properly include "config.h". Perl scripts also now all use #!/usr/bin/env instead of assuming that it's in /usr/bin/perl. (PR#1628, fixes #1606. Reported by Robert Clausecker) * Fixed NAME entry in htslib-s3-plugin man page so the whatis and apropos commands find it. (PR#1634, thanks to Étienne Mollier) * Assorted dependency tracking fixes. (PR#1653, thanks to John Marshall) Documentation updates --------------------- * Changed Alpine build instructions as they've switched back to using openssl. (PR#1609) * Recommend using -rdynamic when statically linking a libhts.a with plugins enabled. (PR#1611, thanks to John Marshall. Fixes #1600, reported by Jack Wimberley) * Fixed example in docs for sam_hdr_add_line(). (PR#1618, thanks to kojix2) * Improved test harness for base modifications API. (PR#1648) Bug fixes --------- * Fix a major bug when searching against a CRAM index where one container has start and end coordinates entirely contained within the previous container. This would occasionally miss data, and sometimes return much more than required. The bug affected versions 1.11 to 1.17, although the change in 1.11 was bug-fixing multi-threaded index queries. This bug did not affect index building. There is no need to reindex your CRAM files. (PR#1574, PR#1640. Fixes #1569, #1639, samtools/samtools#1808, samtools/samtools#1819. Reported by xuxif, Jens Reeder and Jared Simpson) * Prevent CRAM blocks from becoming too big in files with short sequences but very long aux tags. (PR #1613) * Fix bug where the CRAM decoder for CONST_INT and CONST_BYTE codecs may incorrectly look for extra data in the CORE block. Note that this bug only affected the experimental CRAM v4.0 decoder. (PR#1614) * Fix crypt4gh redirection so it works in conjunction with non-file IO, such as using htsget. (PR#1577) * Improve error checking for the VCF POS column, when facing invalid data. (PR#1575, replaces #1570 originally reported and fixed by Colin Nolan.) * Improved error checking on VCF indexing to validate the data is BGZF compressed. (PR#1581) * Fix bug where bin number calculation could overflow when making iterators over regions that go to the end of a chromosome. (PR#1595) * Backport attractivechaos/klib#78 (by Pall Melsted) to HTSlib. Prevents infinite loops in kseq_read() when reading broken gzip files. (PR#1582, fixes #1579. Reported by Goran Vinterhalter) * Backport attractivechaos/klib@384277a (by innoink) to HTSlib. Fixes the kh_int_hash_func2() macro definition. (PR#1599, fixes #1598. Reported by fanxinping) * Remove a compilation warning on systems with newer libcurl releases. (PR#1572) * Windows: Fixed BGZF EOF check for recent MinGW releases. (PR#1601, fixes samtools/bcftools#1901) * Fixed bug where tabix would not return the correct regions for files where the column ordering is end, ..., begin instead of begin, ..., end. (PR#1626, fixes #1622. Reported by Hiruna Samarakoon) * sam_format_aux1() now always NUL-terminates Z/H tags. (PR#1631) * Ensure base modification iterator is reset when no MM tag is present. (PR#1631, PR#1647) * Fix segfault when attempting to write an uncompressed BAM file opened using hts_open(name, "wbu"). This was attempting to write BAM data without wrapping it in BGZF blocks, which is invalid according to the BAM specification. "wbu" is now internally converted to "wb0" to output uncompressed data wrapped in BGZF blocks. (PR#1632, fixes #1617. Reported by Joyjit Daw) * Fixed over-strict bounds check in probaln_glocal() which caused it to make sub-optimal alignments when the requested band width was greater than the query length. (PR#1616, fixes #1605. Reported by Jared Simpson) * Fixed possible double frees when handling errors in bcf_hdr_add_hrec(), if particular memory allocations fail. (PR#1637) * Ensure that bcf_hdr_remove() clears up all pointers to the items removed from dictionaries. Failing to do this could have resulted in a call requesting a deleted item via bcf_hdr_get_hrec() returning a stale pointer. (PR#1637) * Stop the gzip decompresser from finishing prematurely when an empty gzip block is followed by more data. (PR#1643, PR#1646) ------------------------------------------------------------------------------ samtools - changes v1.18 ------------------------------------------------------------------------------ New work and changes: * Add minimiser sort option to collate by an indexed fasta. Expand the minimiser sort to arrange the minimiser values in the same order as they occur in the reference genome. This is acts as an extremely crude and simplistic read aligner that can be used to boost read compression. (PR#1818) * Add a --duplicate-count option to markdup. Adds the number of duplicates (including itself) to the original read in a 'dc' tag. (PR#1816. Thanks to wulj2) * Make calmd handle unaligned data or empty files without throwing an error. This is to make pipelines work more smoothly. A warning will still be issued. (PR#1841, fixes #1839. Reported by Filipe G. Vieira) * Consistent, more comprehensive flag filtering for fasta/fastq. Added --rf/--incl[ude]-flags and long options for -F (--excl[ude]-flags and -f (--require-flags). (PR#1842. Thanks to Devang Thakkar) * Apply fastq --input-fmt-option settings. Previously any options specified were not being applied to the input file. (PR#1855. Thanks to John Marshall) * Add fastq -d TAG[:VAL] check. This mirrors view -d and will only output alignments that match TAG (and VAL if specified). (PR#1863, fixes #1854. Requested by Rasmus Kirkegaard) * Extend import --order TAG to --order TAG:length. If length is specified, the tag format goes from integer to a 0-padded string format. This is a workaround for BAM and CRAM that cannot encode an order tag of over 4 billion records. (PR#1850, fixes #1847. Reported by Feng Tian) * New -aa mode for consensus. This works like the -aa option in depth and mpileup. The single 'a' reports all bases in contigs covered by alignments. Double 'aa' (or '-a -a') reports Ns even for the references with no alignments against them. (PR#1851, fixes #1849. Requested by Tim Fennell) * Add long option support to samtools index. (PR#1872, fixes #1869. Reported by Jason Bacon) * Be consistent with rounding of "average length" in samtools stats. (PR#1876, fixes #1867. Reported by Jelinek-J) * Add option to ampliconclip that marks reads as unmapped when they do not have enough aligned bases left after clipping. Default is to unmap reads with zero aligned bases. (PR#1865, fixes #1856. Requested by ces) Bug Fixes: * [From HTSLib] Fix a major bug when searching against a CRAM index where one container has start and end coordinates entirely contained within the previous container. This would occasionally miss data, and sometimes return much more than required. The bug affected versions 1.11 to 1.17, although the change in 1.11 was bug-fixing multi-threaded index queries. This bug did not affect index building. There is no need to reindex your CRAM files. (PR#samtools/htslib#1574, PR#samtools/htslib#1640. Fixes #samtools/htslib#1569, #samtools/htslib#1639, #1808, #1819. Reported by xuxif, Jens Reeder and Jared Simpson) * Fix a sort -M bug (regression) when merging sub-blocks. Data was valid but in a poor order for compression. (PR#1812) * Fix bug in split output format. Now SAM and CRAM format can chosen as well as BAM. Also a documentation change, see below. (PR#1821) * Add error checking to view -e filter expression code. Invalid expressions were not returning an error code. (PR#1833, fixes #1829. Reported by Steve Huang) * Fix reheader CRAM output version. Sets the correct CRAM output version for non-3.0 CRAMs. (PR#1868, fixes #1866. Reported by John Marshall) Documentation: * Expand the default filtering information on the mpileup man page. (PR#1802, fixes #1801. Reported by gevro) * Add an explanation of the default behaviour of split files on generating a file for reads with missing or unrecognised RG tags. Also a small bug fix, see above. (PR#1821, fixes #1817. Reported by Steve Huang) * In the INSTALL instructions, switched back to openssl for Alpine. This matches the current Alpine Linux practice. (PR#1837, see htslib#1591. Reported by John Marshall) * Fix various typos caught by lintian parsers. (PR#1877. Thanks to Étienne Mollier) * Document consensus --qual-calibration option. (PR#1880, fixes #1879. Reported by John Marshall) * Updated the page about samtools duplicate marking with more detail at www.htslib.org/algorithms/duplicate.html Non user-visible changes and build improvements: * Removed a redundant line that caused a warning in gcc-13. (PR#1838) ------------------------------------------------------------------------------ bcftools - changes v1.18 ------------------------------------------------------------------------------ Changes affecting the whole of bcftools, or multiple commands: * Support auto indexing during writing BCF and VCF.gz via new `--write-index` option Changes affecting specific commands: * bcftools annotate - The `-m, --mark-sites` option can be now used to mark all sites without the need to provide the `-a` file (#1861) - Fix a bug where the `-m` function did not respect the `--min-overlap` option (#1869) - Fix a bug when update of INFO/END results in assertion error (#1957) * bcftools concat - New option `--drop-genotypes` * bcftools consensus - Support higher-ploidy genotypes with `-H, --haplotype` (#1892) - Allow `--mark-ins` and `--mark-snv` with a character, similarly to `--mark-del` * bcftools convert - Support for conversion from tab-delimited files (CHROM,POS,REF,ALT) to sites-only VCFs * bcftools csq - New `--unify-chr-names` option to automatically unify different chromosome naming conventions in the input GFF, fasta and VCF files (e.g. "chrX" vs "X") - More versatility in parsing various flavors of GFF - A new `--dump-gff` option to help with debugging and investigating the internals of hGFF parsing - When printing consequences in nonsense mediated decay transcripts, include 'NMD_transcript' in the consequence part of the annotation. This is to make filtering easier and analogous to VEP annotations. For example the consequence annotation 3_prime_utr|PCGF3|ENST00000430644|NMD is newly printed as 3_prime_utr&NMD_transcript|PCGF3|ENST00000430644|NMD * bcftools gtcheck - Add stats for the number of sites matched in the GT-vs-GT, GT-vs-PL, etc modes. This information is important for interpretation of the discordance score, as only the GT-vs-GT matching can be interpreted as the number of mismatching genotypes. * bcftools +mendelian2 - Fix in command line argument parsing, the `-p` and `-P` options were not functioning (#1906) * bcftools merge - New `-M, --missing-rules` option to control the behavior of merging of vector tags to prevent mixtures of known and missing values in tags when desired - Use values pertaining to the unknown allele (<*> or <NON_REF>) when available to prevent mixtures of known and missing values (#1888) - Revamped line matching code to fix problems in gVCF merging where split gVCF blocks would not update genotypes (#1891, #1164). * bcftool mpileup - Fix a bug in --indels-v2.0 which caused an endless loop when CIGAR operator 'H' or 'P' was encountered * bcftools norm - The `-m, --multiallelics +` mode now preserves phasing (#1893) - Symbolic <DEL.*> alleles are now normalized too (#1919) - New `-g, --gff-annot` option to right-align indels in forward transcripts to follow HGVS 3'rule (#1929) * bcftools query - Force newline character in formatting expression when not given explicitly - Fix `-H` header output in formatting expressions containing newlines * bcftools reheader - Make `-f, --fai` aware of long contigs not representable by 32-bit integer (#1959) * bcftools +split-vep - Prevent a segfault when `-i/-e` use a VEP subfield not included in `-f` or `-c` (#1877) - New `-X, --keep-sites` option complementing the existing `-x, --drop-sites` options - Force newline character in formatting expression when not given explicitly - Fix a subtle ambiguity: identical rows must be returned when `-s` is applied regardless of `-f` containing the `-a` VEP tag itself or not. * bcftools stats - Collect new VAF (variant allele frequency) statistics from FORMAT/AD field - When counting transitions/transversions, consider also alternate het genotypes * plot-vcfstats - Add three new VAF plots -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: James B. <jk...@sa...> - 2023-07-05 15:04:28
|
On Wed, Jul 05, 2023 at 03:11:50PM +0200, Adrian wrote: > I noticed samtools works perfectly fine with gzip SAM and BAM files, but it > does not work (to my knowledge) with gzip CRAM. Is there any reason for > this? Is there any plan to implement support to gzip CRAM directly? CRAM is already superior compression to gzip, so running gzip on the CRAM will give minimal if any benefit. James -- James Bonfield (jk...@sa...) The Sanger Institute, Hinxton, Cambs, CB10 1SA -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is Wellcome Sanger Institute, Wellcome Genome Campus, Hinxton, CB10 1SA. |
From: Adrian <ace...@gm...> - 2023-07-05 13:12:09
|
Hi, I noticed samtools works perfectly fine with gzip SAM and BAM files, but it does not work (to my knowledge) with gzip CRAM. Is there any reason for this? Is there any plan to implement support to gzip CRAM directly? Regards Adrian |
From: Robert D. <rm...@sa...> - 2023-02-21 14:39:24
|
Samtools (and HTSlib and BCFtools) version 1.17 is now available from GitHub and SourceForge. https://github.com/samtools/htslib/releases/tag/1.17 https://github.com/samtools/samtools/releases/tag/1.17 https://github.com/samtools/bcftools/releases/tag/1.17 https://sourceforge.net/projects/samtools/ The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.17 ------------------------------------------------------------------------------ * A new API for iterating through a BAM record's aux field. (PR#1354, addresses #1319. Thanks to John Marshall) * Text mode for bgzip. Allows bgzip to compress lines of text with block breaks at newlines. (PR#1493, thanks to Mike Lin for the initial version PR#1369) * Make tabix support CSI indices with large positions. Unlike SAM and VCF files, BED files do not set a maximum reference length which hindered CSI support. This change sets an arbitrary large size of 100G to enable it to work. (PR#1506) * Add a fai_line_length function. Exposes the internal line-wrap length. (PR#1516) * Check for invalid barcode tags in fastq output. (PR#1518, fixes samtools#1728. Reported by Poshi) * Warn if reference found in a CRAM file is not contained in the specified reference file. (PR#1517 and PR#1521, adds diagnostics for #1515. Reported by Wei WeiDeng) * Add a faidx_seq_len64 function that can return sequence lengths longer than INT_MAX. At the same time limit faidx_seq_len to INT_MAX output. Also add a fai_adjust_region to ensure given ranges do not go beyond the end of the requested sequence. (PR#1519) * Add a bcf_strerror function to give text descriptions of BCF errors. (PR#1510) * Add CRAM SQ/M5 header checking when specifying a fasta file. This is to prevent creating a CRAM that cannot be decoded again. (PR#1522. In response to samtools#1748 though not a direct fix) * Improve support for very long input lines (> 2Gbyte). This is mostly useful for tabix which does not do much interpretation of its input. (PR#1542, a partial fix for #1539) * Speed up load_ref_portion. This function has been sped up by about 7x, which speeds up low-depth CRAM decoding by about 10%. (PR#1551) * Expand CRAM API to cope with new samtools cram_size command. (PR#1546) * Merges neighbouring I and D ops into one op within pileup. This means 4M1D1D1D3M is reported as 4M3D3M. Fixing this in sam.c means not only is samtools mpileup now looking better, but any tool using the mpileup API will be getting consistent results. (PR#1552, fixes the last remaining part of samtools#139) * Update the API documentation for bgzf_mt as it refered to a previous iteration. (PR#1556, fixes #1553. Reported by Raghavendra Padmanabhan) Build changes ------------- * Use POSIX grep in testing as egrep and fgrep are considered obsolete. (PR#1509, thanks to David Seifert) * Switch to building libdefalte with cmake for Cirris CI. (PR#1511) * Ensure strings in config_vars.h are escaped correctly. (PR#1530, fixes #1527. Reported by Lucas Czech) * Easier modification of shared library permissions during install. (PR#1532, fixes #1525. Reported by StephDC) * Fix build on ancient compilers. Added -std=gnu90 to build tests so older C compilers will still be happy. (PR#1524, fixes #1523. Reported by Martin Jakt) * Switch MacOS CI tests to an ARM-based image. (PR#1536) * Cut down the number of embed_ref=2 tests that get run. (PR#1537) * Add symbol versions to libhts.so. This is to aid package developers. (PR#1560 addresses #1505, thanks to John Marshall. Reported by Stefan Bruens) * htscodecs now updated to v1.4.0. (PR#1563) * Cleaned up misleading system error reports in test_bgzf. (PR#1565) Bug fixes --------- * VCF. Fix n-squared complexity in sample line with many adjacent tabs [fuzz]. (PR#1503) * Improved bcftools detection and reporting of bgzf decode errors. (PR#1504, thanks to Lilian Janin. PR#1529 thanks to Bergur Ragnarsson, fixes #1528. PR#1554) * Prevent crash when the only FASTA entry has no sequence [fuzz]. (PR#1507) * Fixed typo in sam.h documentation. (PR#1512, thanks to kojix2) * Fix buffer read-overrun in bam_plp_insertion_mod. (PR#1520) * Fix hash keys being left behind by bcf_hdr_remove. (PR#1535, fixes #1533. Reported by Giulio Genovese in #842) * Make bcf_hdr_idinfo_exists more robust by checking id value exists. (PR#1544, fixes #1538. Reported by Giulio Genovese) * CRAM improvements. Fixed crash with multi-threaded CRAM. Fixed a bug in the codec parameter learning for CRAM 3.1 name tokeniser. Fixed Cram compression container substitution matrix generation, (PR#1558, PR#1559 and PR#1562) ------------------------------------------------------------------------------ samtools - changes v1.17 ------------------------------------------------------------------------------ New work and changes: * New samtools reset subcommand. Removes alignment information. Alignment location, CIGAR, mate mapping and flags are updated. If the alignment was in reverse direction, sequence and its quality values are reversed and complemented and the reverse flag is reset. Supplementary and secondary alignment data are discarded. (PR#1767, implements #1682. Requested by dkj) * New samtools cram-size subcommand. It writes out metrics about a CRAM file reporting aggregate sizes per block "Content ID" fields, the data-series contained within them, and the compression methods used. (PR#1777) * Added a --sanitize option to fixmate and view. This performs some sanity checks on the state of SAM record fields, fixing up common mistakes made by aligners. (PR#1698) * Permit 1 thread with samtools view. All other subcommands already allow this and it does provide a modest speed increase. (PR#1755, fixes #1743. Reported by Goran Vinterhalter) * Add CRAM_OPT_REQUIRED_FIELDS option for view -c. This is a big speed up for CRAM (maybe 5-fold), but it depends on which filtering options are being used. (PR#1776, fixes #1775. Reported by Chang Y) * New filtering options in samtools depth. The new --excl-flags option is a synonym for -G, with --incl-flags and --require-flags added to match view logic. (PR#1718, fixes #1702. Reported by Dario Beraldi) * Speed up calmd's slow handling of non-position-sorted data by adding caching. This uses more memory but is only activated when needed. (PR#1723, fixes #1595. Reported by lxwgcool) * Improve samtools consensus for platforms with instrument specific profiles, considerably helping for data with very different indel error models and providing base quality recalibration tables. On PacBio HiFi, ONT and Ultima Genomics consensus qualities are also redistributed within homopolymers and the likelihood of nearby indel errors is raised. (PR#1721, PR#1733) * Consensus --mark-ins option. This permits he consensus output to include a markup indicating the next base is an insertion. This is necessary as we need a way of outputting both consensus and also how that consensus marries up with the reference coordinates. (PR#1746) * Make faidx/fqidx output line length default to the input line length. (PR#1738, fixes #1734. Reported by John Marshall) * Speed up optical duplicate checking where data has a lot of duplicates compared to non-duplicates. (PR#1779, fixes #1771. Reported by Poshi) * For collate use TMPDIR environment variable, when looking for a temporary folder. (PR#1782, based on PR#1178 and fixes #1172. Reported by Martin Pollard) Bug Fixes: * Fix stats breakage on long deletions when given a reference. (PR#1712, fixes #1707. Reported by John Didion) * In ampliconclip, stop hard clipping from wrongly removing entire reads. (PR#1722, fixes #1717. Reported by Kevin Xu) * Fix bug in ampliconstats where references mentioned in the input file headers but not in the bed file would cause it to complain that the SAM headers were inconsistent. (PR#1727, fixes #1650. Reported by jPontix) * Fixed SEGV in samtools collate when no filename given. (PR#1724) * Changed the default UMI barcode regex in markdup. The old regex was too restrictive. This version will at least allow the default read name UMI as given in the Illumina example documentation. (PR#1737, fixes #1730. Reported by yloemie) * Fix samtools consensus buffer overrun with MD:Z handling. (PR#1745, fixes #1744. Reported by trilisser) * Fix a buffer read-overflow in mpileup and tview on sequences with seq "*". (PR#1747) * Fix view -X command line parsing that was broken in 1.15. (PR#1772, fixes #1720. Reported by Francisco Rodríguez-Algarra and Miguel Machado) * Stop samtools view -d from reporting meaningless system errors when tag validation fails. (PR#1796) Documentation: * Add a description of the samtools tview display layout to the man page. Documents . vs , and upper vs lowercase. Adds a -s sample example, and documents the -w option. (PR#1765, fixes #1759. Reported by Lucas Ferreira da Silva) * Clarify intention of samtools fasta/q in man page and soft vs hard clipping. (PR#1794, fixes #1792. Reported by Ryan Lorig-Roach) * Minor fix to wording of mpileup --rf usage and man page. (PR#1795, fixes #1791. Reported by Luka Pavageau) Non user-visible changes and build improvements: * Use POSIX grep in testing as egrep and fgrep are considered obsolete. (PR#1726, thanks to David Seifert) * Switch MacOS CI tests to an ARM-based image. (PR#1770) ------------------------------------------------------------------------------ bcftools - changes v1.17 ------------------------------------------------------------------------------ Changes affecting the whole of bcftools, or multiple commands: * The -i/-e filtering expressions - Error checks were added to prevent incorrect use of vector arithmetics. For example, when evaluating the sum of two vectors A and B, the resulting vector could contain nonsense values when the input vectors were not of the same length. The fix introduces the following logic: - evaluate to C_i = A_i + B_i when length(A)==B(A) and set length(C)=length(A) - evaluate to C_i = A_i + B_0 when length(B)=1 and set length(C)=length(A) - evaluate to C_i = A_0 + B_i when length(A)=1 and set length(C)=length(B) - throw an error when length(A)!=length(B) AND length(A)!=1 AND length(B)!=1 - Arrays in Number=R tags can be now subscripted by alleles found in FORMAT/GT. For example, FORMAT/AD[GT] > 10 .. require support of more than 10 reads for each allele FORMAT/AD[0:GT] > 10 .. same as above, but in the first sample sSUM(FORMAT/AD[GT]) > 20 .. require total sample depth bigger than 20 * The commands `consensus -H` and `+split-vep -H` - Drop unnecessary leading space in the first header column and newly print `#[1]columnName` instead of the previous `# [1]columnName` (#1856) Changes affecting specific commands: * bcftools +allele-length - Fix overflow for indels longer than 512bp and aggregate alleles equal or larger than that in the same bin (#1837) * bcftools annotate - Support sample reordering of annotation file (#1785) - Restore lost functionality of the --pair-logic option (#1808) * bcftools call - Fix a bug where too many alleles passed to `-C alleles` via `-T` caused memory corruption (#1790) - Fix a bug where indels constrained with `-C alleles -T` would sometimes be missed (#1706) * bcftools consensus - BREAKING CHANGE: the option `-I, --iupac-codes` newly outputs IUPAC codes based on FORMAT/GT of all samples. The `-s, --samples` and `-S, --samples-file` options can be used to subset samples. In order to ignore samples and consider only the REF and ALT columns (the original behavior prior to 1.17), run with `-s -` (#1828) * bcftools convert - Make variantkey conversion work for sites without an ALT allele (#1806) * bcftool csq - Fix a bug where a MNV with multiple consequences (e.g. missense + stop_gained) would report only the less severe one (#1810) - GFF file parsing was made slightly more flexible, newly ids can be just 'XXX' rather than, for example, 'gene:XXX' - New gff2gff perl script to fix GFF formatting differences * bcftools +fill-tags - More of the available annotations are now added by the `-t all` option * bcftools +fixref - New INFO/FIXREF annotation - New -m swap mode * bcftools +mendelian - The +mendelian plugin has been deprecated and replaced with +mendelian2. The function of the plugin is the same but the command line options and the output format has changed, and for this was introduced as a new plugin. * bcftools mpileup - Most of the annotations generated by mpileup are now optional via the `-a, --annotate` option and add several new (mostly experimental) annotations. - New option `--indels-2.0` for an EXPERIMENTAL indel calling model. This model aims to address some known deficiencies of the current indel calling algorithm, specifically, it uses diploid reference consensus sequence. Note that in the current version it has the potential to increase sensitivity but at the cost of decreased specificity. - Make the FS annotation (Fisher exact test strand bias) functional and remove it from the default annotations * bcftools norm - New --multi-overlaps option allows to set overlapping alleles either to the ref allele (the current default) or to a missing allele (#1764 and #1802) - Fixed a bug in `-m -` which does not split missing FORMAT values correctly and could lead to empty FORMAT fields such as `::` instead of the correct `:.:` (#1818) - The `--atomize` option previously would not split complex indels such as C>GGG. Newly these will be split into two records C>G and C>CGG (#1832) * bcftools query - Fix a rare bug where the printing of SAMPLE field with `query` was incorrectly suppressed when the `-e` option contained a sample expression while the formatting query did not. See #1783 for details. * bcftools +setGT - Add new `--new-gt X` option (#1800) - Add new `--target-gt r:FLOAT` option to randomly select a proportion of genotypes (#1850) - Fix a bug where `-t ./x` mode was advertised as selecting both phased and unphased half-missing genotypes, but was in fact selecting only unphased genotypes (#1844) * bcftools +split-vep - New options `-g, --gene-list` and `--gene-list-fields` which allow to prioritize consequences from a list of genes, or restrict output to the listed genes - New `-H, --print-header` option to print the header with `-f` - Work around a bug in the LOFTEE VEP plugin used to annotate gnomAD VCFs. There the LoF_info subfield contains commas which, in general, makes it impossible to parse the VEP subfields. The +split-vep plugin can now work with such files, replacing the offending commas with slash (/) characters. See also https://github.com/Ensembl/ensembl-vep/issues/1351 - Newly the `-c, --columns` option can be omitted when a subfield is used in `-i/-e` filtering expression. Note that `-c` may still have to be given when it is not possible to infer the type of the subfield. Note that this is an experimental feature. * bcftools stats - The per-sample stats (PSC) would not be computed when `-i/-e` filtering options and the `-s -` option were given but the expression did not include sample columns (1835) * bcftools +tag2tag - Revamp of the plugin to allow wider range of tag conversions, specifically all combinations from FORMAT/GL,PL,GP to FORMAT/GL,PL,GP,GT * bcftools +trio-dnm2 - New `-n, --strictly-novel` option to downplay alleles which violate Mendelian inheritance but are not novel - Allow to set the `--pn` and `--pns` options separately for SNVs and indels and make the indel settings more strict by default - Output missing FORMAT/VAF values in non-trio samples, rather than random nonsense values * bcftools +variant-distance - New option `-d, --direction` to choose the directionality: forward, reverse, nearest (the default) or both (#1829) -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: gilsanan <gil...@i3...> - 2022-11-29 19:05:29
|
Hi, I'm new using samtools and I'm sorry in advance if my question has an obvious answer... I did not manage to find any solution online or in the samtools documentation. I've been trying to filter the reads in a BAM file by the reference name. This field appears in the 3rd column of the header in the bam file. I checked my BAM file and it seems correct. This is an example of one read in my file: > f5f22aad-bb4a-4ec6-805c-4f290cf5d116 16 MA1#0001 51 60 35S188M1D11M76S * 0 0 ATTGCTAATGATACGGCGACCACCGAGATCTACACGTTCAGAGTTCTACAGTCCGACGATCCATTGTTCACCTCACGCAGACGAGCTCTTAGACAAGTGACAGTGTCCACTCGAGACTGGATATTCTCCAGTTACATACCAATAATAACAAATACTTAACATTGTACAAGAATACCCGATTCTGCGTCCAACCTCATTGAAACCCTCCTTTGGGTGGAATTCTGGGTGCCAAGGAATCCAGTCACCACGATATCTCGTATGCCGTCTTCTGCTTGAGCAATACGTAACTGAACGAAGTGCCGACAAAAGA 4>?ABEEEDC@@@?@>>>>>A@?ABBCAAAABCBCCDBBBDBDCBABA@@ABCCCC?>>>>HDDCDAAAA@@@@@@BCBBAABCCD?==<<?@ABCBC@@??>?@BBICBBABFHFDDDIQHHEEDCCDCBBCBCDDDEGFHHGEDDEGIEEEEEECCC@@AACLEFCBBBBAAABCGIKB@?>>?B@>>>>ADBCCCA?;;;<>ACC?>>?A?>??<=,++(78;:2-,+(''((+-2CBAB?>???FCB?>?@?DEC::::D///,-...-?>GEFDAAAAH6555;50,*'%$#"""##$$####"" NM:i:1 ms:i:392 AS:i:392 nn:i:0 tp:A:P cm:i:31 s1:i:168 s2:i:51 de:f:0.005 rl:i:18 where MA1#0001 is the name of one of my references. I was trying to filter the file by reference name with no luck (always empty outputs). I checked the statistics then using the following line: > samtools idxstats my_file.bam The output shows 0 reads mapped to every reference, here a fragment of the output: > MA1#0001 300 1099 0 > MA1#0002 300 2325 0 > MA1#0003 300 3757 0 > MA1#0004 300 1686 0 > MA1#0005 300 1856 0 > MA1#0006 300 3051 0 > MA1#0007 300 2965 0 Any idea about why this is happening an how to solve it? Many thanks!! Cheers, Eva |
From: Trine R. <tri...@uk...> - 2022-11-09 12:29:04
|
Hi, I have without much luck tried to install samtools on my pc. I am using WSL2. when I type the command "make" during installation I get the following error: /usr/bin/ld: cannot open output file samtools: Is a directory collect2: error: ld returned 1 exit status make: *** [Makefile:145: samtools] Error 1 Does someone have a solution as to what my problem could be?, Best, Trine |
From: James B. <jk...@sa...> - 2022-10-25 14:56:30
|
On Tue, Oct 25, 2022 at 04:25:32PM +0200, Carla Mariner Llicer wrote: > I'm using Samtools version 1.3.1. I was analysing some illumina paired-end Please use a modern version. 1.3 is 6 year old, with many bug fixes and improvements made since then. This particular issue was fixed in 1.13. However, also see wish-list item https://github.com/samtools/samtools/issues/1354 for potential issues still associated with random selection vs by-quality. James -- James Bonfield (jk...@sa...) The Sanger Institute, Hinxton, Cambs, CB10 1SA -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: Carla M. L. <cma...@gm...> - 2022-10-25 14:25:51
|
Hi! I'm using Samtools version 1.3.1. I was analysing some illumina paired-end sequenced samples and I realised that when running mpileup with this version, by default, samtools discards all the Minus Strand reads when there are overlapping reads. Is there any reason why always discards Minus Strand reads? Does it make sense? What happens is that I have some real SNPs in regions where all the reads overlap so, when it discards all the Minus Strand reads, the SNPs appear only in one strand and then it doesn't pass the downstream filters that we apply. I know that with the option --ignore-overlaps I can keep all the overlapping reads so my question is Why samtools does not discard Positive strands or Minus Strands in a random way by default? Thanks for your attention, Carla |
From: James B. <jk...@sa...> - 2022-09-08 08:08:19
|
On Wed, Sep 07, 2022 at 09:27:17PM +0000, Leal, Joel via Samtools-help wrote: > Thanks for your answer. Does it mean I need to install v1.16.1 instead of v1.2? Yes. See https://github.com/samtools/samtools/releases for the list of what's changed with each release. A quick search for consensus in that page shows it appeared with 1.15. Wherever you read that it was in 1.2 was simply wrong. James -- James Bonfield (jk...@sa...) The Sanger Institute, Hinxton, Cambs, CB10 1SA -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: Tom B. <tb...@um...> - 2022-09-07 21:28:01
|
Joel - I strongly suspect that the problem is "Which is the higher version number, samtools v1.16.1 or samtools v1.2 ?" In this case it is v1.16.1. Our sysadmins have had trouble with this too and have "updated" us back to version 1.2. - thanks - tom blackwell - On Wed, Sep 7, 2022 at 4:40 PM Leal, Joel via Samtools-help < sam...@li...> wrote: > Hi, > > I am trying to generate a consensus sequence from a .BAM file using > Samtools. I found it was developed for Samtools1.2. However, when calling > the function consensus, it is not recognized. The function mpileup is > recognized. > > Could you tell me what is the issue, please? > > Thanks, > Joel > > > _______________________________________________ > Samtools-help mailing list > Sam...@li... > https://lists.sourceforge.net/lists/listinfo/samtools-help > |
From: Leal, J. <jd...@he...> - 2022-09-07 21:27:35
|
Thanks for your answer. Does it mean I need to install v1.16.1 instead of v1.2? Joel ________________________________ From: Tom Blackwell <tb...@um...> Sent: Wednesday, September 7, 2022 2:01 PM To: Leal, Joel <jd...@he...> Cc: sam...@li... <sam...@li...> Subject: Re: [Samtools-help] consensus not recognized Joel - I strongly suspect that the problem is "Which is the higher version number, samtools v1.16.1 or samtools v1.2 ?" In this case it is v1.16.1. Our sysadmins have had trouble with this too and have "updated" us back to version 1.2. - thanks - tom blackwell - On Wed, Sep 7, 2022 at 4:40 PM Leal, Joel via Samtools-help <sam...@li...<mailto:sam...@li...>> wrote: Hi, I am trying to generate a consensus sequence from a .BAM file using Samtools. I found it was developed for Samtools1.2. However, when calling the function consensus, it is not recognized. The function mpileup is recognized. Could you tell me what is the issue, please? Thanks, Joel [cid:18319bef83bcb971f161] _______________________________________________ Samtools-help mailing list Sam...@li...<mailto:Sam...@li...> https://lists.sourceforge.net/lists/listinfo/samtools-help<https://urldefense.com/v3/__https://lists.sourceforge.net/lists/listinfo/samtools-help__;!!LLK065n_VXAQ!kj4Uu5XaZVIF-MsGG3-Os_mgJQ7N-eLd0ETcyVrgDMB-rnlkpvf829NXdb_16i7peqfKgOr5GhwbhkHrd_g$> |
From: Leal, J. <jd...@he...> - 2022-09-07 20:40:11
|
Hi, I am trying to generate a consensus sequence from a .BAM file using Samtools. I found it was developed for Samtools1.2. However, when calling the function consensus, it is not recognized. The function mpileup is recognized. Could you tell me what is the issue, please? Thanks, Joel [cid:ac963de6-540b-40c8-a8c8-0cfd7d336023] |
From: Robert D. <rm...@sa...> - 2022-09-02 14:20:17
|
Samtools version 1.16.1 is now available from GitHub and SourceForge. This release fixes some bugs in the new template-coordinate sort feature. https://sourceforge.net/projects/samtools/ https://github.com/samtools/samtools/releases/tag/1.16.1 The main changes are listed below: ------------------------------------------------------------------------------ samtools - changes v1.16.1 ------------------------------------------------------------------------------ Bug fixes: * Fixed a bug with the template-coordinate sort which caused incorrect ordering when using threads, or processing large files that don't fit completely in memory. (PR#1703, thanks to Nils Homer) * Fixed a crash that occurred when trying to use `samtools merge` in template-coordinate mode. (PR#1705, thanks to Nils Homer) -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: Robert D. <rm...@sa...> - 2022-08-18 14:17:46
|
Samtools (and HTSlib and BCFtools) version 1.16 is now available from GitHub and SourceForge. https://sourceforge.net/projects/samtools/ https://github.com/samtools/htslib/releases/tag/1.16 https://github.com/samtools/samtools/releases/tag/1.16 https://github.com/samtools/bcftools/releases/tag/1.16 The main changes are listed below: ------------------------------------------------------------------------------ htslib - changes v1.16 ------------------------------------------------------------------------------ * Make hfile_s3 refresh AWS credentials on expiry in order to make HTSlib work better with AWS IAM credentials, which have a limited lifespan. (PR#1462 and PR#1474, addresses #344) * Allow BAM headers between 2GB and 4GB in size once more. This is not permitted in the BAM specification but was allowed in an earlier version of HTSlib. There is now a warning at 2GB and a hard failure at 4GB. (PR#1421, fixes #1420 and samtools#1613. Reported by John Marshall and R C Mueller) * Improve error message when failing to load an index. (PR#1468, example of the problem samtools#1637) * Permit MM (base modification) tags containing "." and "?" suffixes. These define implicit vs explicit coordinates. See the SAM tags specification for details. (PR#1423 and PR#1426, fixes #1418. PR#1469, fixes #1466. Reported by cjw85) * Warn if spaces instead of tabs are detected in a VCF file to prevent confusion. (PR#1328, fixes bcftools#1575. Reported by ketkijoshi278) * Add an "sclen" filter expression keyword. This is the length of a soft-clip, both left and right end. It may be combined with qlen (qlen-sclen) to obtain the number of bases in the query sequence that have been aligned to the genome ie it provides a way to compare local-alignment vs global-alignment length. (PR#1441 and PR/samtools#1661, fixes #1436. Requested by Chang Y) * Improve error messages for CRAM reference mismatches. If the user specifies the wrong reference, the CRAM slice header MD5sum checks fail. We now report the SQ line M5 string too so it is possible to validate against the whole chr in the ref.fa file. The error message has also been improved to report the reference name instead of #num. Finally, we now hint at the likely cause, which counters the misleading samtools supplied error of "truncated or corrupt" file. (PR#1427, fixes samtools#1640. Reported by Jian-Guo Zhou) * Expose more of the CRAM API and add new functionality to extract the reference from a CRAM file. (PR#1429 and PR#1442) * Improvements to the implementation of embedded references in CRAM where no external reference is specified. (PR#1449, addresses some of the issues in #1445) * The CRAM writer now allows alignment records with RG:Z: aux tags that don't have a corresponding @RG ID in the file header. Previously these tags would have been silently dropped. HTSlib will complain whenever it has to add one though, as such tags do not conform to recommended practice for the SAM, BAM and CRAM formats. (PR#1480, fixes #1479. Reported by Alex Leonard) * Set tab delimiter in man page for tabix GFF3 sort. (PR#1457. Thanks to Colin Diesh) * When using libdeflate, the 1...9 scale of BGZF compression levels is now remapped to the 1...12 range used by libdeflate instead of being passed directly. In particular, HTSlib levels 8 and 9 now map to libdeflate levels 10 and 12, so it is possible to select the highest (but slowest) compression offered by libdeflate. (PR#1488, fixes #1477. Reported by Gert Hulselmans) * The VCF variant API has been extended so that it can return separate flags for INS and DEL variants as well as the existing INDEL one. These flags have not been added to the old bcf_get_variant_types() interface as it could break existing users. To access them, it is necessary to use new functions bcf_has_variant_type() and bcf_has_variant_types(). (PR#1467) * The missing, but trivial, `le_to_u8()` function has been added to hts_endian. (PR#1494, Thanks to John Marshall) * bcf_format_gt() now works properly on big-endian platforms. (PR#1495, Thanks to John Marshall) Build changes ------------- These are compiler, configuration and makefile based changes. * Update htscodecs to version 1.3.0 for new SIMD code + various fixes. Updates the htscodecs submodule and adds changes necessary to make HTSlib build the new SIMD codec implementations. (PR#1438, PR#1489, PR#1500) * Fix clang builds under mingw. Under mingw, clang requires dllexport to be applied to both function declarations and function definitions. (PR#1435, PR#1497, PR#1498 fixes #1433. Reported by teepean) * Fix curl type warning with gcc 12.1 on Windows. (PR#1443) * Detect ARM Neon support and only build appropriate SIMD object files. (PR#1451, fixes #1450. Thanks to John Marshall) * `make print-config` now reports extra CFLAGS that are needed to build the SIMD parts of htscodecs. These may be of use to third-party build systems that don't use HTSlib's or htscodecs' build infrastructure. (PR#1485. Thanks to John Marshall) * Fixed some Makefile dependency issues for the "check"/"test" targets and plugins. In particular, "make check" will now build the "all" target, if not done already, before running the tests. (PR#1496) Bug fixes --------- * Fix bug when reading position -1 in BCF (0 in VCF), which is used to indicate telomeric regions. The BCF reader was incorrectly assuming the value stored in the file was unsigned, so a VCF->BCF->VCF round-trip would change it from 0 to 4294967296. (PR#1476, fixes #1475 and bcftools#1753. Reported by Rodrigo Martin) * Various bugs and quirks have been fixed in the filter expression engine, mostly related to the handling of absent tags, and the is_true flag. Note that as a result of these fixes, some filter expressions may give different results: - Fixed and-expressions including aux tag values which could give an invalid true result depending on the order of terms. - The expression `![NM]` is now true if only `NM` does not exist. In earlier versions it would also report true for tags like `NM:i:0` which exist but have a value of zero. - The expression `[X1] != 0` is now false when `X1` does not exist. Earlier versions would return true for this comparison when the tag was missing. - NULL values due to missing tags now propagate through string, bitwise and mathematical operations. Logical operations always treat them as false. (PR#1463, fixes samtools#1670. Reported by Gert Hulselmans; PR#1478, fixes samtools#1677. Reported by johnsonzcode) * Fix buffer overrun in bam_plp_insertion_mod. Memory now grows to the proper size needed for base modification data. (PR#1430, fixes samtools#1652. Reported by hd2326) * Remove limit of returned size from fai_retrieve(). (PR#1446, fixes samtools#1660. Reported by Shane McCarthy) * Cap hts_getline() return value at INT_MAX. Prevents hts_getline() from returning a negative number (a fail) for very long string length values. (PR#1448. Thanks to John Marshall) * Fix breakend detection and test bcf_set_variant_type(). (PR#1456, fixes #1455. Thanks to Martin Pollard) * Prevent arrays of BCF_BT_NULL values found in BCF files from causing bcf_fmt_array() to call exit() as the type is unsupported. These are now tested for and caught by bcf_record_check(), which returns an error code instead. (PR#1486) * Improved detection of fasta and fastq files that have very long comments following identifiers. (PR#1491, thanks to John Marshall. Fixes samtools/samtools#1689, reported by cjw85) * Fixed a SEGV triggered by giving a SAM file to `samtools import`. (PR#1492) ------------------------------------------------------------------------------ samtools - changes v1.16 ------------------------------------------------------------------------------ New work and changes: * samtools reference command added. This subcommand extracts the embedded reference out of a CRAM file. (PR#1649, addresses #723. Requested by Torsten Seemann) * samtools import now adds grouped by query-name to the header. (PR#1633, thanks to Nils Homer) * Made samtools view read error messages more generic. Former error message would claim that there was a "truncated file or corrupt BAM index file" with no real justification. Also reset errno in stream_view which could lead to confusing error messages. (PR#1645, addresses some of the issues in #1640. Reported by Jian-Guo Zhou) * Make samtools view -p also clear mqual, tlen and cigar. (PR#1647, fixes #1606. Reported by eboyden) * Add bedcov option -c to report read count. (PR#1644, fixes #1629. Reported by Natchaphon Rajudom) * Add UMI/barcode handling to samtools markdup. (PR#1630, fixes #1358 and #1514. Reported by Gert Hulselmans and Poshi) * Add a new template coordinate sort order to samtools sort and samtools merge. This is useful when working with unique molecular identifiers (UMIs). (PR#1605, fixes #1591. Thanks to Nils Homer) * Rename mpileup --ignore-overlaps to --ignore-overlaps-removal or --disable-overlap-removal. The previous name was ambiguous and was often read as an option to enable removal of overlapping bases, while in reality this is on by default and the option turns off the ability to remove overlapping bases. (PR#1666, fixes #1663. Reported by yangdingyangding) * The dict command can now read BWA's .alt file and add AH:* tags indicating reference sequences that represent alternate loci. (PR#1676. Thanks to John Marshall) * The "samtools index" command can now accept multiple alignment filenames with the new -M option, and will index each of them separately. (Specifying the output index filename via out.index or the new -o option is currently only applicable when there is only one alignment file to be indexed.) (PR#1674. Reported by Abigail Ramsøe and Nicola Romanò. Thanks to John Marshall) * Allow samtools fastq -T "*". This allows all tags from SAM records to be written to fastq headers. This is a counterpart to samtools import -T "*". (PR#1679. Thanks to cjw85) Bug Fixes: * Re-enable --reference option for samtools depth. The reference is not used but this makes the command line usage compatible with older releases. (PR#1646, fixes #1643. Reported by Randy Harr) * Fix regex coordinate bug in samtools markdup. (PR#1657, fixes #1642. Reported by Randy Harr) * Fix divide by zero in plot-bamstats -m, on unmapped data. (PR#1678, fixes #1675. Thanks to Shane McCarthy) * Fix missing RG headers when using samtools merge -r. (PR#1683, addresses htslib#1479. Reported by Alex Leonard) * Fix a possible unaligned access in samtools reference. (PR#1696) Documentation: * Add documentation on CRAM compression profiles and some of the newer options that appear in CRAM 3.1 and above. (PR#1659, fixes #1656. Reported by Matthias De Smet) * Add "sclen" filter expression keyword documentation. (PR#1661, see also htslib#1441) * Extend FILTER EXPRESSION man page section to match the changes made in HTSlib. (PR#1687, samtools/htslib#1478) Non user-visible changes and build improvements: * Ensure generated test files are ignored (by git) and cleaned (by make testclean) (PR#1692, Thanks to John Marshall) ------------------------------------------------------------------------------ bcftools - changes v1.16 ------------------------------------------------------------------------------ * New plugin `bcftools +variant-distance` to annotate records with distance to the nearest variant (#1690) Changes affecting the whole of bcftools, or multiple commands: * The -i/-e filtering expressions - Added support for querying of multiple filters, for example `-i 'FILTER="A;B"'` can be used to select sites with two filters "A" and "B" set. See the documentation for more examples. - Added modulo arithmetic operator Changes affecting specific commands: * bcftools annotate - A bug introduced in 1.14 caused that records with INFO/END annotation would incorrectly trigger `-c ~INFO/END` mode of comparison even when not explicitly requested, which would result in not transferring the annotation from a tab-delimited file (#1733) * bcftools merge - New `-m snp-ins-del` switch to merge SNVs, insertions and deletions separately (#1704) * bcftools mpileup - New NMBZ annotation for Mann-Whitney U-z test on number of mismatches within supporting reads - Suppress the output of MQSBZ and FS annotations in absence of alternate allele * bcftools +scatter - Fix erroneous addition of duplicate PG lines * bcftools +setGT - Custom genotypes (e.g. `-n c:1/1`) now correctly override ploidy -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: James B. <jk...@sa...> - 2022-08-04 08:19:28
|
On Thu, Aug 04, 2022 at 08:43:58AM +0100, Thomas Juettemann wrote: > Hi Rob, > Thanks for looking into it. Unfortunately isec keeps only the first record. If you're not after a full intersection of only-in-A only-in-B and in-both, then it's possible you could use filtering options instead. Eg "bcftools view -T A.vcf.gz B.vcf.gz" will report records from B that overlap locations listed in A. It doesn't need to a be BED file as it'll auto-detect the file format. > > On Tue, 2 Aug 2022, Thomas Juettemann wrote: > > > > > I came across a "transcript-based" VCF file, meaning a variant can be > > > present multiple times but belonging to a different transcript. See > > > "FIle 1" below as an example. I am finding myself in the unfortunate > > > situation of having to intersect ("File 2") and retain all records > > > with the same position and REF/ALT ("Desired output"). > > > Long shot: Is that possible? > > > > Does "bcftools isec" (https://www.htslib.org/doc/bcftools.html#isec) do > > what you want? The "Extract and write records from A shared by both A and > > B using exact allele match" example in the manual page sounds like it > > might: > > > > bcftools isec -p dir -n=2 -w1 A.vcf.gz B.vcf.gz I think this means that in the above example, multiple transcripts in B that overlap the coordinates in A will still be shown. If you need the reverse, then it'd need another command with A and B swapped around. I'm not sure this is exactly the same thing, but it's worth an experiment with a few simple examples to validate it. (Take care with complex variants and not just SNPs to check "overlap" works as you expect when indels are present.) James -- James Bonfield (jk...@sa...) The Sanger Institute, Hinxton, Cambs, CB10 1SA -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |
From: Thomas J. <jue...@gm...> - 2022-08-04 07:44:21
|
Hi Rob, Thanks for looking into it. Unfortunately isec keeps only the first record. Digging through Github, it seems it is a known limitation: https://github.com/samtools/bcftools/issues/665#issuecomment-323372893 Best, Thomas On Wed, 3 Aug 2022 at 17:50, Robert Davies <rm...@sa...> wrote: > > On Tue, 2 Aug 2022, Thomas Juettemann wrote: > > > I came across a "transcript-based" VCF file, meaning a variant can be > > present multiple times but belonging to a different transcript. See > > "FIle 1" below as an example. I am finding myself in the unfortunate > > situation of having to intersect ("File 2") and retain all records > > with the same position and REF/ALT ("Desired output"). > > Long shot: Is that possible? > > Does "bcftools isec" (https://www.htslib.org/doc/bcftools.html#isec) do > what you want? The "Extract and write records from A shared by both A and > B using exact allele match" example in the manual page sounds like it > might: > > bcftools isec -p dir -n=2 -w1 A.vcf.gz B.vcf.gz > > If not, you can't find anything else, and you only want to do a few of > them, it might be possible to break out pysam and write something. If you > want to do lots, then a C program would probably be the way forward - it > doesn't look like it would be too difficult. > > Rob Davies rm...@sa... > The Sanger Institute http://www.sanger.ac.uk/ > Hinxton, Cambs., Tel. +44 (1223) 834244 > CB10 1SA, U.K. Fax. +44 (1223) 494919 > > > -- > The Wellcome Sanger Institute is operated by Genome Research > Limited, a charity registered in England with number 1021457 and a > company registered in England with number 2742969, whose registered > office is 215 Euston Road, London, NW1 2BE. |
From: Robert D. <rm...@sa...> - 2022-08-03 16:51:06
|
On Tue, 2 Aug 2022, Thomas Juettemann wrote: > I came across a "transcript-based" VCF file, meaning a variant can be > present multiple times but belonging to a different transcript. See > "FIle 1" below as an example. I am finding myself in the unfortunate > situation of having to intersect ("File 2") and retain all records > with the same position and REF/ALT ("Desired output"). > Long shot: Is that possible? Does "bcftools isec" (https://www.htslib.org/doc/bcftools.html#isec) do what you want? The "Extract and write records from A shared by both A and B using exact allele match" example in the manual page sounds like it might: bcftools isec -p dir -n=2 -w1 A.vcf.gz B.vcf.gz If not, you can't find anything else, and you only want to do a few of them, it might be possible to break out pysam and write something. If you want to do lots, then a C program would probably be the way forward - it doesn't look like it would be too difficult. Rob Davies rm...@sa... The Sanger Institute http://www.sanger.ac.uk/ Hinxton, Cambs., Tel. +44 (1223) 834244 CB10 1SA, U.K. Fax. +44 (1223) 494919 -- The Wellcome Sanger Institute is operated by Genome Research Limited, a charity registered in England with number 1021457 and a company registered in England with number 2742969, whose registered office is 215 Euston Road, London, NW1 2BE. |