Answers.Body: <p>I think what you are looking for is to perform a <a href="https://www.nature.com/articles/s43586-021-00056-9" rel="nofollow noreferrer">GWAS</a>: Genome-Wide Association Study. Generally, these are used to detect links between mutations in patients and complex diseases. I see no reason that the same approach couldn't be used to attempt to predict someone's height.</p>
  <p>I haven't run a GWAS in a long time, so I wouldn't be able to recommend any specific tools to use, but I found packages in <a href="https://cran.r-project.org/web/packages/statgenGWAS/vignettes/GWAS.html" rel="nofollow noreferrer">R</a> and in <a href="https://hail.is/docs/0.2/tutorials/01-genome-wide-association-study.html" rel="nofollow noreferrer">Python</a> that could work.</p>
  <p>One thing you didn't mention is how you know the height of each of your donors, which would seem like the biggest complexity to me. Once you have the group VCF + associated NGS files, you're mostly covered from a data standpoint.</p>
  <p>If you're looking to go deeper, you could investigate dimension reduction methods like PCA that can also serve to associate variables with certain traits, or even some machine learning methods like SVMs or a Neural-Net regression model.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2024-06-18T01:46:53.627Z
Answers.Text: Thanks for the response! My target is actually plant, it was the plant height I have been referring to, but I guess that is not making much differences? There were several samples involved in the sampling and genotyping. The height of plants were measured at nth years. I am confused with what you have mentioned: "group VCF + associated NGS files", could you elaborate more on this?
Answers.UserDisplayName: web
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2024-06-18T08:48:57.483Z
Answers.Text: GWAS should work on plants just fine, even if it was designed for use on Human diseases (these two articles use GWAS on plants https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8316211/ https://www.nature.com/articles/s41396-021-00993-z maybe their methods detail their parametres). The "group VCF" refers to the file containing the variant calls for your entire population, which can be quite long to obtain and requires a lot of computing power, so having access to it is a boon. If you also have access to the raw NGS data, then you can also explore other methods of analysis, like k-mers, if you want.
Answers.UserDisplayName: Whitehot
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2024-06-13T15:09:18.073Z
Answers.LastActivityDate: 2024-06-14T13:00:02.650Z
Answers.LastEditDate: 2024-06-14T13:00:02.650Z
Answers.UserDisplayName: Whitehot
Body: <p>I'm new to handling SNP genotyping data in VCF format.</p>
  <p>My goal is to identify significant SNPs linked to a specific trait (like &quot;height&quot;) for specific samples in a multi-sample VCF file. My previous colleague who initially handled this project left, and I'm struggling to replicate their workflow for filtering and trait-based SNP extraction.</p>
  <p>Here's what I'm working with: a large VCF file containing data from 32 samples, ideally needing a filtered subset focused on just 5 specific samples. These VCF files are GATK outputs with all detected SNPs. I want to extract statistically significant SNPs that are non-homogeneous (meaning they vary) for the target trait (height).</p>
  <p>GATK command in the VCF file (If this helps):</p>
  <pre><code>##GATKCommandLine=&lt;ID=SelectVariants,Version=3.3-0-g37228af,Date=&quot;Thu Jan 12 15:24:24 CET 2017&quot;,Epoch=1484231064776,CommandLineOptions=&quot;analysis_type=SelectVariants input_file=[] showFullBamList=false read_buffer_size=null phone_home=AWS gatk_key=null tag=NA read_filter=[] intervals=null excludeIntervals=null interval_set_rule=UNION interval_merging=ALL interval_padding=0 reference_sequence=/work/NGSwaiting4newNAS/G01_chromosome.fna nonDeterministicRandomSeed=false disableDithering=false maxRuntime=-1 maxRuntimeUnits=MINUTES downsampling_type=BY_SAMPLE downsample_to_fraction=null downsample_to_coverage=1000 baq=OFF baqGapOpenPenalty=40.0 refactor_NDN_cigar_string=false fix_misencoded_quality_scores=false allow_potentially_misencoded_quality_scores=false useOriginalQualities=false defaultBaseQualities=-1 performanceLog=null BQSR=null quantize_quals=0 disable_indel_quals=false emit_original_quals=false preserve_qscores_less_than=6 globalQScorePrior=-1.0 validation_strictness=SILENT remove_program_records=false keep_program_records=false sample_rename_mapping_file=null unsafe=null disable_auto_index_creation_and_locking_when_reading_rods=false no_cmdline_in_header=false sites_only=false never_trim_vcf_format_field=false bcf=false bam_compression=null simplifyBAM=false disable_bam_indexing=false generate_md5=false num_threads=1 num_cpu_threads_per_data_thread=1 num_io_threads=0 monitorThreadEfficiency=false num_bam_file_handles=null read_group_black_list=null pedigree=[] pedigreeString=[] pedigreeValidationType=STRICT allow_intervals_with_unindexed_bam=false generateShadowBCF=false variant_index_type=DYNAMIC_SEEK variant_index_parameter=-1 logging_level=INFO log_to_file=null help=false version=false variant=(RodBinding name=variant source=DBG_SNP_EGRESQ_filtered.vcf) discordance=(RodBinding name= source=UNBOUND) concordance=(RodBinding name= source=UNBOUND) out=org.broadinstitute.gatk.engine.io.stubs.VariantContextWriterStub sample_name=[] sample_expressions=null sample_file=null exclude_sample_name=[] exclude_sample_file=[] select_expressions=[] excludeNonVariants=false excludeFiltered=false restrictAllelesTo=BIALLELIC keepOriginalAC=false mendelianViolation=false mendelianViolationQualThreshold=0.0 select_random_fraction=0.0 remove_fraction_genotypes=0.0 selectTypeToInclude=[SNP] keepIDs=null fullyDecode=false forceGenotypesDecode=false justRead=false maxIndelSize=2147483647 ALLOW_NONOVERLAPPING_COMMAND_LINE_SAMPLES=false filter_reads_with_N_cigar=false filter_mismatching_base_and_quals=false filter_bases_not_stored=false&quot;&gt;
  </code></pre>
  <p>Further complicating things, I found previously filtered VCFs with significant SNPs for 5 samples. Now, I'm tasked with repeating the process for more samples.</p>
  <p>While I've researched VCF filtering by origin selection, I'm stuck on pinpointing the most efficient workflow to achieve both filtering and trait-based SNP extraction. Any recommendations for tools or resources would be incredibly helpful! Please feel free to ask if there is any elaboration needed. Thank you!</p>
  <p><strong>Update:</strong> I also have another file with a list of QTL regions that were significant to the height trait.</p>
  <p>Below is a picture generated to help with understanding of the challenges:
  <a href="https://i.sstatic.net/Z98kK9mS.jpg" rel="nofollow noreferrer"><img src="https://i.sstatic.net/Z98kK9mS.jpg" alt="Figure" /></a>
  I have information on the specific region of QTL that are significant to the height trait. But the information is based on all the available origins. Now I am only interested to some of the origins. How do I filter the VCF files to get the SNP related to the interested origins?</p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-06T12:39:30.320Z
Comments.Text: What does "origin" mean? Do you mean different individuals were sequenced and all results are in a single, multi-sample VCF? And do you _only_ have SNPs? Nothing larger? Do you have a budget or do you need to only use free tools?
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-07T03:09:50.640Z
Comments.Text: @terdon Yes, origin refers to different sample included in the genotyping and all the results were compiled in one VCF. I also have the QTL regions that are significant to the trait, the fasta sequences for each chromosomes, etc. I think the files are complete as they have completed the task previously. I would like to explore the free tools if possible.
Comments.UserDisplayName: web
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-07T10:37:16.880Z
Comments.Text: OK, and what do you have in the VCF file? Is that annotated with the gene the variant falls on/near? Do you have a list of genes for your phenotype? Again, are these only SNPs (normal GATK output would also include indels)?
Comments.UserDisplayName: terdon
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-07T13:00:25.553Z
Comments.Text: `extract statistically significant SNPs that are non-homogeneous (meaning they vary) for the target trait (height).`  assuming controls (A, B, C) and samples with the desired trait (X, Y, Z) what does it mean? In the order of the letters above, anything different from say 6 Ts, ABC=reference base, XYZ=T, ABC=whatever, XYZ=CGT?
Comments.UserDisplayName: darked89
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-10T01:48:13.547Z
Comments.Text: @terdon The SNP is not annotated with the gene, I think the information is available somewhere else, but not in the VCF file. At the moment, I can only see list of QTL significant to the phenotype and several sequences related to the phenotype from publication. I am not sure how to check if indels were included in the VCF file, is it like what explained by Martijin here (https://www.biostars.org/p/48731/#48736)? If yes, then no, as the ref and alt allele consisted only a single base pair here?
Comments.UserDisplayName: web
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-10T01:56:41.637Z
Comments.Text: @darked89 From what I understand, in the VCF file, there are a reference and alternate allele. The genotyping result for each sample are based on this. For example the genotype received by sample A is 0/0, sample B is 0/1 and sample C is 1/1. In this case, sample A is having G/G, sample C is having T/T, which are homogeneous, only sample B is having G/T, the non-homogeneous I mean in the text.
Comments.UserDisplayName: web
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-10T02:24:45.073Z
Comments.Text: @terdon are the indels usually being noted in the "ID" column? My "ID" columns are filled with only dot. I have updated the question with GATK command I found in the VCF file.
Comments.UserDisplayName: web
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2024-06-10T13:13:18.017Z
Comments.Text: No, the ID column will be empty unless you've annotated the VCF with a tool that populates them. Indels appear as you said: >1 residue in either the REF or ALT fields, or both.
Comments.UserDisplayName: terdon
Tags: <vcf><genome><snp><variant-calling><gatk>
Title: Replicating VCF Filtering & Trait-Based SNP Extraction Workflow
ThreadText: # Replicating VCF Filtering & Trait-Based SNP Extraction Workflow
  
  (Asked by: username_0 on 2024-06-06)
  
  I'm new to handling SNP genotyping data in VCF format.
  
  My goal is to identify significant SNPs linked to a specific trait (like "height") for specific samples in a multi-sample VCF file. My previous colleague who initially handled this project left, and I'm struggling to replicate their workflow for filtering and trait-based SNP extraction.
  
  Here's what I'm working with: a large VCF file containing data from 32 samples, ideally needing a filtered subset focused on just 5 specific samples. These VCF files are GATK outputs with all detected SNPs. I want to extract statistically significant SNPs that are non-homogeneous (meaning they vary) for the target trait (height).
  
  GATK command in the VCF file (If this helps):
  
  ```
  ##GATKCommandLine=<ID=SelectVariants,Version=3.3-0-g37228af,Date="Thu Jan 12 15:24:24 CET 2017",Epoch=1484231064776,CommandLineOptions="analysis_type=SelectVariants input_file=[] showFullBamList=false read_buffer_size=null phone_home=AWS gatk_key=null tag=NA read_filter=[] intervals=null excludeIntervals=null interval_set_rule=UNION interval_merging=ALL interval_padding=0 reference_sequence=/work/NGSwaiting4newNAS/G01_chromosome.fna nonDeterministicRandomSeed=false disableDithering=false maxRuntime=-1 maxRuntimeUnits=MINUTES downsampling_type=BY_SAMPLE downsample_to_fraction=null downsample_to_coverage=1000 baq=OFF baqGapOpenPenalty=40.0 refactor_NDN_cigar_string=false fix_misencoded_quality_scores=false allow_potentially_misencoded_quality_scores=false useOriginalQualities=false defaultBaseQualities=-1 performanceLog=null BQSR=null quantize_quals=0 disable_indel_quals=false emit_original_quals=false preserve_qscores_less_than=6 globalQScorePrior=-1.0 validation_strictness=SILENT remove_program_records=false keep_program_records=false sample_rename_mapping_file=null unsafe=null disable_auto_index_creation_and_locking_when_reading_rods=false no_cmdline_in_header=false sites_only=false never_trim_vcf_format_field=false bcf=false bam_compression=null simplifyBAM=false disable_bam_indexing=false generate_md5=false num_threads=1 num_cpu_threads_per_data_thread=1 num_io_threads=0 monitorThreadEfficiency=false num_bam_file_handles=null read_group_black_list=null pedigree=[] pedigreeString=[] pedigreeValidationType=STRICT allow_intervals_with_unindexed_bam=false generateShadowBCF=false variant_index_type=DYNAMIC_SEEK variant_index_parameter=-1 logging_level=INFO log_to_file=null help=false version=false variant=(RodBinding name=variant source=DBG_SNP_EGRESQ_filtered.vcf) discordance=(RodBinding name= source=UNBOUND) concordance=(RodBinding name= source=UNBOUND) out=org.broadinstitute.gatk.engine.io.stubs.VariantContextWriterStub sample_name=[] sample_expressions=null sample_file=null exclude_sample_name=[] exclude_sample_file=[] select_expressions=[] excludeNonVariants=false excludeFiltered=false restrictAllelesTo=BIALLELIC keepOriginalAC=false mendelianViolation=false mendelianViolationQualThreshold=0.0 select_random_fraction=0.0 remove_fraction_genotypes=0.0 selectTypeToInclude=[SNP] keepIDs=null fullyDecode=false forceGenotypesDecode=false justRead=false maxIndelSize=2147483647 ALLOW_NONOVERLAPPING_COMMAND_LINE_SAMPLES=false filter_reads_with_N_cigar=false filter_mismatching_base_and_quals=false filter_bases_not_stored=false">
  
  ```
  
  Further complicating things, I found previously filtered VCFs with significant SNPs for 5 samples. Now, I'm tasked with repeating the process for more samples.
  
  While I've researched VCF filtering by origin selection, I'm stuck on pinpointing the most efficient workflow to achieve both filtering and trait-based SNP extraction. Any recommendations for tools or resources would be incredibly helpful! Please feel free to ask if there is any elaboration needed. Thank you!
  
  **Update:** I also have another file with a list of QTL regions that were significant to the height trait.
  
  Below is a picture generated to help with understanding of the challenges:
  [![Figure](https://i.sstatic.net/Z98kK9mS.jpg)](https://i.sstatic.net/Z98kK9mS.jpg)
  I have information on the specific region of QTL that are significant to the height trait. But the information is based on all the available origins. Now I am only interested to some of the origins. How do I filter the VCF files to get the SNP related to the interested origins?
  
  **Tags:** vcf, genome, snp, variant-calling, gatk
  
  > **Comments:**
  > **username_1**: What does "origin" mean? Do you mean different individuals were sequenced and all results are in a single, multi-sample VCF? And do you _only_ have SNPs? Nothing larger? Do you have a budget or do you need to only use free tools?  
  > **username_0**: @username_1 Yes, origin refers to different sample included in the genotyping and all the results were compiled in one VCF. I also have the QTL regions that are significant to the trait, the fasta sequences for each chromosomes, etc. I think the files are complete as they have completed the task previously. I would like to explore the free tools if possible.  
  > **username_1**: OK, and what do you have in the VCF file? Is that annotated with the gene the variant falls on/near? Do you have a list of genes for your phenotype? Again, are these only SNPs (normal GATK output would also include indels)?  
  > **username_2**: `extract statistically significant SNPs that are non-homogeneous (meaning they vary) for the target trait (height).`  assuming controls (A, B, C) and samples with the desired trait (X, Y, Z) what does it mean? In the order of the letters above, anything different from say 6 Ts, ABC=reference base, XYZ=T, ABC=whatever, XYZ=CGT?  
  > **username_0**: @username_1 The SNP is not annotated with the gene, I think the information is available somewhere else, but not in the VCF file. At the moment, I can only see list of QTL significant to the phenotype and several sequences related to the phenotype from publication. I am not sure how to check if indels were included in the VCF file, is it like what explained by Martijin here (https://www.biostars.org/p/48731/#48736)? If yes, then no, as the ref and alt allele consisted only a single base pair here?  
  > **username_0**: @username_2 From what I understand, in the VCF file, there are a reference and alternate allele. The genotyping result for each sample are based on this. For example the genotype received by sample A is 0/0, sample B is 0/1 and sample C is 1/1. In this case, sample A is having G/G, sample C is having T/T, which are homogeneous, only sample B is having G/T, the non-homogeneous I mean in the text.  
  > **username_0**: @username_1 are the indels usually being noted in the "ID" column? My "ID" columns are filled with only dot. I have updated the question with GATK command I found in the VCF file.  
  > **username_1**: No, the ID column will be empty unless you've annotated the VCF with a tool that populates them. Indels appear as you said: >1 residue in either the REF or ALT fields, or both.  
  
  ---
  
  ## 1 answers
  
  ### Answer by username_3 on 2024-06-13. Score: 2
  
  I think what you are looking for is to perform a [GWAS](https://www.nature.com/articles/s43586-021-00056-9): Genome-Wide Association Study. Generally, these are used to detect links between mutations in patients and complex diseases. I see no reason that the same approach couldn't be used to attempt to predict someone's height.
  
  I haven't run a GWAS in a long time, so I wouldn't be able to recommend any specific tools to use, but I found packages in [R](https://cran.r-project.org/web/packages/statgenGWAS/vignettes/GWAS.html) and in [Python](https://hail.is/docs/0.2/tutorials/01-genome-wide-association-study.html) that could work.
  
  One thing you didn't mention is how you know the height of each of your donors, which would seem like the biggest complexity to me. Once you have the group VCF + associated NGS files, you're mostly covered from a data standpoint.
  
  If you're looking to go deeper, you could investigate dimension reduction methods like PCA that can also serve to associate variables with certain traits, or even some machine learning methods like SVMs or a Neural-Net regression model.
  
  > **Comments:**
  > **username_0**: Thanks for the response! My target is actually plant, it was the plant height I have been referring to, but I guess that is not making much differences? There were several samples involved in the sampling and genotyping. The height of plants were measured at nth years. I am confused with what you have mentioned: "group VCF + associated NGS files", could you elaborate more on this?  
  > **username_3**: GWAS should work on plants just fine, even if it was designed for use on Human diseases (these two articles use GWAS on plants https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8316211/ https://www.nature.com/articles/s41396-021-00993-z maybe their methods detail their parametres). The "group VCF" refers to the file containing the variant calls for your entire population, which can be quite long to obtain and requires a lot of computing power, so having access to it is a boon. If you also have access to the raw NGS data, then you can also explore other methods of analysis, like k-mers, if you want.  
  
  ---
  
  