{"action":"created","issue":{"url":"https://api.github.com/repos/broadinstitute/gatk/issues/3487","repository_url":"https://api.github.com/repos/broadinstitute/gatk","labels_url":"https://api.github.com/repos/broadinstitute/gatk/issues/3487/labels{/name}","comments_url":"https://api.github.com/repos/broadinstitute/gatk/issues/3487/comments","events_url":"https://api.github.com/repos/broadinstitute/gatk/issues/3487/events","html_url":"https://github.com/broadinstitute/gatk/issues/3487","id":252308599,"node_id":"MDU6SXNzdWUyNTIzMDg1OTk=","number":3487,"title":"Port LeftAlignAndTrimVariants and improve stdout messaging...","user":{"login":"sooheelee","id":11543866,"node_id":"MDQ6VXNlcjExNTQzODY2","avatar_url":"https://avatars2.githubusercontent.com/u/11543866?v=4","gravatar_id":"","url":"https://api.github.com/users/sooheelee","html_url":"https://github.com/sooheelee","followers_url":"https://api.github.com/users/sooheelee/followers","following_url":"https://api.github.com/users/sooheelee/following{/other_user}","gists_url":"https://api.github.com/users/sooheelee/gists{/gist_id}","starred_url":"https://api.github.com/users/sooheelee/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/sooheelee/subscriptions","organizations_url":"https://api.github.com/users/sooheelee/orgs","repos_url":"https://api.github.com/users/sooheelee/repos","events_url":"https://api.github.com/users/sooheelee/events{/privacy}","received_events_url":"https://api.github.com/users/sooheelee/received_events","type":"User","site_admin":false},"labels":[],"state":"closed","locked":false,"assignee":{"login":"sooheelee","id":11543866,"node_id":"MDQ6VXNlcjExNTQzODY2","avatar_url":"https://avatars2.githubusercontent.com/u/11543866?v=4","gravatar_id":"","url":"https://api.github.com/users/sooheelee","html_url":"https://github.com/sooheelee","followers_url":"https://api.github.com/users/sooheelee/followers","following_url":"https://api.github.com/users/sooheelee/following{/other_user}","gists_url":"https://api.github.com/users/sooheelee/gists{/gist_id}","starred_url":"https://api.github.com/users/sooheelee/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/sooheelee/subscriptions","organizations_url":"https://api.github.com/users/sooheelee/orgs","repos_url":"https://api.github.com/users/sooheelee/repos","events_url":"https://api.github.com/users/sooheelee/events{/privacy}","received_events_url":"https://api.github.com/users/sooheelee/received_events","type":"User","site_admin":false},"assignees":[{"login":"sooheelee","id":11543866,"node_id":"MDQ6VXNlcjExNTQzODY2","avatar_url":"https://avatars2.githubusercontent.com/u/11543866?v=4","gravatar_id":"","url":"https://api.github.com/users/sooheelee","html_url":"https://github.com/sooheelee","followers_url":"https://api.github.com/users/sooheelee/followers","following_url":"https://api.github.com/users/sooheelee/following{/other_user}","gists_url":"https://api.github.com/users/sooheelee/gists{/gist_id}","starred_url":"https://api.github.com/users/sooheelee/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/sooheelee/subscriptions","organizations_url":"https://api.github.com/users/sooheelee/orgs","repos_url":"https://api.github.com/users/sooheelee/repos","events_url":"https://api.github.com/users/sooheelee/events{/privacy}","received_events_url":"https://api.github.com/users/sooheelee/received_events","type":"User","site_admin":false},{"login":"kachulis","id":39926576,"node_id":"MDQ6VXNlcjM5OTI2NTc2","avatar_url":"https://avatars1.githubusercontent.com/u/39926576?v=4","gravatar_id":"","url":"https://api.github.com/users/kachulis","html_url":"https://github.com/kachulis","followers_url":"https://api.github.com/users/kachulis/followers","following_url":"https://api.github.com/users/kachulis/following{/other_user}","gists_url":"https://api.github.com/users/kachulis/gists{/gist_id}","starred_url":"https://api.github.com/users/kachulis/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/kachulis/subscriptions","organizations_url":"https://api.github.com/users/kachulis/orgs","repos_url":"https://api.github.com/users/kachulis/repos","events_url":"https://api.github.com/users/kachulis/events{/privacy}","received_events_url":"https://api.github.com/users/kachulis/received_events","type":"User","site_admin":false}],"milestone":null,"comments":9,"created_at":"2017-08-23T14:43:14Z","updated_at":"2018-09-06T19:56:25Z","closed_at":"2018-09-06T19:56:25Z","author_association":"CONTRIBUTOR","body":"...because it does the following:\r\n- makes it appear no variants were left-aligned when in fact many were left-aligned. \r\n- In addition, because the tool emits an INFO message for sites that it skips (due to reference allele being too long), it appears doubly confusing as if the tool skips sites that could use left-alignment and only these types of skipped sites are present.\r\n- Finally, for my file, it drops 66,568 variant sites (I did not check alleles) for me silently and without any warning.\r\n\r\n## Two suggestions for improvement for LeftAlign's stdout:\r\n1. A count of variants/sites that were left-aligned OR at least NOT the confusing `0 variants were aligned` message.\r\n2. The length of the longest reference allele that the tool skipped.\r\n3. At least the number if not the sites of dropped variant records.\r\n\r\nMy workaround in the meanwhile: I had to `md5` headerless VCFs to confirm the input is different from the output of the tool. \r\n\r\n---\r\n## Here's how I came to the above\r\nQuick note because I use GATK3.7's LeftAlignAndTrimVariants below. I am using LeftAlignAndTrimVariants from GATK3.7.\r\n\r\n### Command\r\nI am using a grch37-->grch38 liftedover (via Picard) VCF as input. It would be more ideal if Picard LiftoverVcf would left-align variants but we cannot have everything we want when we want it. So here I am asking for a better LeftAlign.\r\n\r\nThe output file is for use with GATK4-Mutect2.\r\n\r\n```\r\njava -jar $GATK -T LeftAlignAndTrimVariants -R ~/Documents/ref/hg38/Homo_sapiens_assembly38.fasta -V zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz -o eta_af-only-gnomad_Hg19toGRCh38.vcf.gz\r\n```\r\n\r\n### Portions of stdout\r\n```\r\nINFO  17:39:56,364 HelpFormatter - -------------------------------------------------------------------------------- \r\nINFO  17:39:56,365 HelpFormatter - The Genome Analysis Toolkit (GATK) v3.7-0-gcfedb67, Compiled 2016/12/12 11:21:18 \r\nINFO  17:39:56,366 HelpFormatter - Copyright (c) 2010-2016 The Broad Institute \r\nINFO  17:39:56,366 HelpFormatter - For support and documentation go to https://software.broadinstitute.org/gatk \r\nINFO  17:39:56,366 HelpFormatter - [Tue Aug 22 17:39:56 EDT 2017] Executing on Mac OS X 10.11.6 x86_64 \r\nINFO  17:39:56,366 HelpFormatter - Java HotSpot(TM) 64-Bit Server VM 1.8.0_111-b14 \r\nINFO  17:39:56,368 HelpFormatter - Program Args: -T LeftAlignAndTrimVariants -R /Users/shlee/Documents/ref/hg38/Homo_sapiens_assembly38.fasta -V zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz -o eta_af-only-gnomad_Hg19toGRCh38.vcf.gz \r\nINFO  17:39:56,373 HelpFormatter - Executing as shlee@WMCF9-CB5 on Mac OS X 10.11.6 x86_64; Java HotSpot(TM) 64-Bit Server VM 1.8.0_111-b14. \r\nINFO  17:39:56,374 HelpFormatter - Date/Time: 2017/08/22 17:39:56 \r\nINFO  17:39:56,374 HelpFormatter - -------------------------------------------------------------------------------- \r\nINFO  17:39:56,374 HelpFormatter - -------------------------------------------------------------------------------- \r\nINFO  17:39:56,385 GenomeAnalysisEngine - Strictness is SILENT \r\nINFO  17:39:57,377 GenomeAnalysisEngine - Downsampling Settings: Method: BY_SAMPLE, Target Coverage: 1000 \r\nWARN  17:39:57,688 IndexDictionaryUtils - Track variant doesn't have a sequence dictionary built in, skipping dictionary validation \r\nINFO  17:39:58,497 GenomeAnalysisEngine - Preparing for traversal \r\nINFO  17:39:58,502 GenomeAnalysisEngine - Done preparing for traversal \r\nINFO  17:39:58,503 ProgressMeter - [INITIALIZATION COMPLETE; STARTING PROCESSING] \r\nINFO  17:39:58,503 ProgressMeter -                 | processed |    time |    per 1M |           |   total | remaining \r\nINFO  17:39:58,503 ProgressMeter -        Location |     sites | elapsed |     sites | completed | runtime |   runtime \r\nINFO  17:39:58,692 LeftAlignAndTrimVariants - Reference allele is too long (245) at position chr1:10146; skipping that record. Set --reference_window_stop >= 245 \r\nINFO  17:39:58,697 LeftAlignAndTrimVariants - Reference allele is too long (225) at position chr1:10178; skipping that record. Set --reference_window_stop >= 225 \r\nINFO  17:39:58,700 LeftAlignAndTrimVariants - Reference allele is too long (221) at position chr1:10213; skipping that record. Set --reference_window_stop >= 221 \r\nINFO  17:39:58,700 LeftAlignAndTrimVariants - Reference allele is too long (223) at position chr1:10218; skipping that record. Set --reference_window_stop >= 223 \r\n...\r\nINFO  21:11:24,340 LeftAlignAndTrimVariants - Reference allele is too long (275) at position chr1:213151784; skipping that record. Set --reference_window_stop >= 275 \r\nINFO  21:11:24,341 LeftAlignAndTrimVariants - Reference allele is too long (226) at position chr1:213151832; skipping that record. Set --reference_window_stop >= 226 \r\nINFO  21:11:24,489 LeftAlignAndTrimVariants - Reference allele is too long (202) at position chr1:213428401; skipping that record. Set --reference_window_stop >= 202 \r\nINFO  21:11:25,060 LeftAlignAndTrimVariants - Reference allele is too long (242) at position chr1:214470660; skipping that record. Set --reference_window_stop >= 242 \r\nINFO  21:11:25,060 LeftAlignAndTrimVariants - Reference allele is too long (242) at position chr1:214470687; skipping that record. Set --reference_window_stop >= 242 \r\nINFO  21:11:25,060 LeftAlignAndTrimVariants - Reference allele is too long (242) at position chr1:214470688; skipping that record. Set --reference_window_stop >= 242 \r\nINFO  21:11:25,060 LeftAlignAndTrimVariants - Reference allele is too long (242) at position chr1:214470696; skipping that record. Set --reference_window_stop >= 242 \r\nINFO  21:11:25,061 LeftAlignAndTrimVariants - Reference allele is too long (242) at position chr1:214470746; skipping that record. Set --reference_window_stop >= 242 \r\nINFO  21:11:25,195 LeftAlignAndTrimVariants - Reference allele is too long (260) at position chr1:214679630; skipping that record. Set --reference_window_stop >= 260 \r\nINFO  21:11:25,196 LeftAlignAndTrimVariants - Reference allele is too long (260) at position chr1:214679633; skipping that record. Set --reference_window_stop >= 260 \r\nINFO  21:11:25,196 LeftAlignAndTrimVariants - Reference allele is too long (254) at position chr1:214679636; skipping that record. Set --reference_window_stop >= 254 \r\nINFO  21:11:25,196 LeftAlignAndTrimVariants - Reference allele is too long (248) at position chr1:214679642; skipping that record. Set --reference_window_stop >= 248 \r\nINFO  21:11:25,196 LeftAlignAndTrimVariants - Reference allele is too long (236) at position chr1:214679654; skipping that record. Set --reference_window_stop >= 236 \r\n...\r\nINFO  21:38:54,227 LeftAlignAndTrimVariants - Reference allele is too long (212) at position chr2_KI270894v1_alt:202602; skipping that record. Set --reference_window_stop >= 212 \r\nINFO  21:38:54,233 LeftAlignAndTrimVariants - Reference allele is too long (220) at position chr2_KI270894v1_alt:204859; skipping that record. Set --reference_window_stop >= 220 \r\nINFO  21:38:54,237 LeftAlignAndTrimVariants - Reference allele is too long (262) at position chr2_KI270894v1_alt:207863; skipping that record. Set --reference_window_stop >= 262 \r\n0 variants were aligned\r\nINFO  21:38:54,554 ProgressMeter -            done   3.31246907E8    31.8 m       5.0 s       99.7%    31.8 m       5.0 s \r\nINFO  21:38:54,554 ProgressMeter - Total runtime 1905.29 secs, 31.75 min, 0.53 hours \r\n------------------------------------------------------------------------------------------\r\nDone. There were 4 WARN messages, the first 4 are repeated below.\r\nWARN  17:39:57,688 IndexDictionaryUtils - Track variant doesn't have a sequence dictionary built in, skipping dictionary validation \r\nWARN  18:13:42,039 SimpleTimer - Clock drift of -1,503,348,737,016,211,299 - -1,503,346,772,578,127,937 = 1,964,438,083,362 nanoseconds detected, vs. max allowable drift of 5,000,000,000. Assuming checkpoint/restart event. \r\nWARN  20:14:18,043 SimpleTimer - Clock drift of -1,503,355,916,564,964,097 - -1,503,348,737,015,111,124 = 7,179,549,852,973 nanoseconds detected, vs. max allowable drift of 5,000,000,000. Assuming checkpoint/restart event. \r\nWARN  21:10:35,064 SimpleTimer - Clock drift of -1,503,359,203,412,549,926 - -1,503,355,916,564,817,209 = 3,286,847,732,717 nanoseconds detected, vs. max allowable drift of 5,000,000,000. Assuming checkpoint/restart event. \r\n------------------------------------------------------------------------------------------\r\nWMCF9-CB5:Mutect2 shlee$ \r\n```\r\n\r\n### Notice the following line from above\r\n\r\n> 0 variants were aligned\r\n\r\nAlso, it would be great if the tool, which appears to keep track of the lengths of reference alleles that are too long, could give me the **maximum length** reference allele so that I can go back and set the `--reference_window_stop` argument appropriately in a second round so that I can left-align _all_ of my variants. \r\n\r\n### MD5 and looking into the files, we see input and output are different and in fact the tool did change allele representations:\r\n```\r\nWMCF9-CB5:Mutect2 shlee$ gzcat zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz | grep -v '##' > zeta_headless.txt\r\nWMCF9-CB5:Mutect2 shlee$ md5 zeta_headless.txt \r\nMD5 (zeta_headless.txt) = 6d93f1ea32c99ae5020881fa8265bdc1\r\n\r\nWMCF9-CB5:Mutect2 shlee$ gzcat eta_af-only-gnomad_Hg19toGRCh38.vcf.gz | grep -v '##' > eta_headless.txt\r\nWMCF9-CB5:Mutect2 shlee$ md5 eta_headless.txt \r\nMD5 (eta_headless.txt) = 1316ffa278316f1850e3fd30979d0422\r\n```\r\n\r\nWhen I peak into the file I see left-alignment.\r\n**zeta**\r\n```\r\nchr1    10143   .       CTAACCCCT       C,CTTAACCCCT    6393.44 PASS    AC=12,4;AF=0.000305,0.0001017\r\nchr1    10144   .       TAACCCCTAACCCTAACCCTAACCC       CAACCCCTAACCCTAACCCTAACCC,TACCCCTAACCCTAACCCTAACCC,TTAACCCTAACCCTAACCC  82960.9 PASS    AC=13,16,9;AF=0.0003364,0.0004141,0.0002329\r\nchr1    10145   .       AACCCCT TACCCCT,ACCCT   51129.9 PASS    AC=6,17;AF=0.0001555,0.0004407\r\n```\r\n**eta**\r\n```\r\nchr1    10143   .       CTAACCCCT       C,CTTAACCCCT    6393.44 PASS    AC=12,4;AF=0.000305,0.0001017\r\nchr1    10144   .       TAACCCC CAACCCC,TACCCC,T        82960.90        PASS    AC=13,16,9;AF=0.0003364,0.0004141,0.0002329\r\nchr1    10145   .       AAC     TAC,A   51129.90        PASS    AC=6,17;AF=0.0001555,0.0004407\r\n```\r\n\r\n---\r\nI spot checked three of the 66568 sites that were dropped silently. They are ugly. Given these three are in regions of repeats, perhaps they are not vital to include in the scheme of somatic calling. However, I would have liked to know that there were a number of sites dropped during the leftalign for my accounting.\r\n\r\n```\r\nWMCF9-CB5:Mutect2 shlee$ gzcat zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz | grep '985368\\tPASS\\tAC=5380,12;AF=0.141,0.0003145'\r\nchr1\t10146\t.\t\r\nACCCCTAACCCTAACCCTAACCCTAACCCTAACCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC\t\r\nACCCTAACCCTAACCCTAACCCTAACCCTAACCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC,ACCTAACCCTAACCCTAACCCTAACCCTAACCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC\t985368\tPASS\tAC=5380,12;AF=0.141,0.0003145\r\n\r\nWMCF9-CB5:Mutect2 shlee$ gzcat zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz | grep '2651.95\\tPASS\\tAC=2;AF=5.307e-05'\r\nchr1\t10213\t.\t\r\nCCCTAACCCTAACCCTAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAA\t\r\nCCCCTAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAA\t2651.95\tPASS\tAC=2;AF=5.307e-05\r\n\r\nWMCF9-CB5:Mutect2 shlee$ gzcat zeta_af-only-gnomad_Hg19toGRCh38.vcf.gz | grep '58703.1\\tPASS\\tAC=7;AF=0.0002078'\r\nchr1\t10229\t.\t\r\nAACCCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC\t\r\nACCCTAACCCTAACCCTAAACCCTAAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCAACCCCAACCCCAACCCCAACCCCAACCCCAACCCTAACCCCTAACCCTAACCCTAACCCTACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCCTAACCCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC\t58703.1\tPASS\tAC=7;AF=0.0002078\r\n```\r\n\r\n---\r\nIn conclusion, I could do without the confusing stdout messaging. In fact, we can do better and provide more informative messaging seemingly without any cost to compute."},"comment":{"url":"https://api.github.com/repos/broadinstitute/gatk/issues/comments/419221218","html_url":"https://github.com/broadinstitute/gatk/issues/3487#issuecomment-419221218","issue_url":"https://api.github.com/repos/broadinstitute/gatk/issues/3487","id":419221218,"node_id":"MDEyOklzc3VlQ29tbWVudDQxOTIyMTIxOA==","user":{"login":"sooheelee","id":11543866,"node_id":"MDQ6VXNlcjExNTQzODY2","avatar_url":"https://avatars2.githubusercontent.com/u/11543866?v=4","gravatar_id":"","url":"https://api.github.com/users/sooheelee","html_url":"https://github.com/sooheelee","followers_url":"https://api.github.com/users/sooheelee/followers","following_url":"https://api.github.com/users/sooheelee/following{/other_user}","gists_url":"https://api.github.com/users/sooheelee/gists{/gist_id}","starred_url":"https://api.github.com/users/sooheelee/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/sooheelee/subscriptions","organizations_url":"https://api.github.com/users/sooheelee/orgs","repos_url":"https://api.github.com/users/sooheelee/repos","events_url":"https://api.github.com/users/sooheelee/events{/privacy}","received_events_url":"https://api.github.com/users/sooheelee/received_events","type":"User","site_admin":false},"created_at":"2018-09-06T19:56:25Z","updated_at":"2018-09-06T19:56:25Z","author_association":"CONTRIBUTOR","body":"I am happy with the port and the improvements to the tool by @kachulis. Many thanks. "}}
{"actor":{"display_login":"sooheelee"}}