Answers.Body: <p>The genomes under analysis are under predominant linkage disequilibrium, i.e. linkage. Chromosome 1 shows strong (or stronger - cut-off depending) HW signal as do one or two telomeres - telomeres are known to undergo recombination, so this may not be unsurprising. Nothing wrong with the HW analysis for the core conclusions, its the stuff between the red and blue lines thats difficult (Q-Q plots).</p>
  <p>The Q-Q plots are a bit more complex and I would go with the lower graph you presented together with its HW-analysis. There is clearly part of the data you are not accounting for, @gringer says the Q-Q plot conform to the beta distribution .... they should be bang linear. Lets look in more detail,</p>
  <ul>
  <li>0 to 2 ... predominant disequilibrium that looks fine, nice and linear</li>
  <li>5.5+ .... HW signal that looks ok, nice and linear</li>
  <li>2 to 5.5 ... complicated</li>
  </ul>
  <p>Given you have a cut-off (red line) at around 5.5 then I'd sort of try and talk your way out of this one. You can do a full investigation, such as transformation to linearise the entire Q-Q plot. Personally I've never done this within HW analysis so I dunno, I've done it for GLM. I'd just say they between the blue and red line represents a complex HW signal which we can't be interpreted but either side of that boundary we're good to go. In other words keep it simple.</p>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-20T23:59:35.253Z
Answers.Text: "gringer says the Q-Q plot conform to the beta distribution" - no, I didn't say this. The *p-values* should conform to the linear distribution if the association is not omnigenic, but the *test statistic* [which is not a p-value] will not fit a linear distribution. Most GWAS programs produce beta values as their association statistic, which is why I mentioned beta.
Answers.UserDisplayName: gringer
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T06:53:10.387Z
Answers.Text: Thanks Michael for your opinion. I am coming from the ML field were evaluation of a model with training data is more clear to me. I found it very hard to estimate GWAS filtered by "feeling" and looking at plots (also because I want to do it automatically). I am new to GWAS and all I need is a simple GWAS. I was unsure about using PCA or LMM  or I could go for a simple linear analysis and whether this looks kind of "normal".Thanks for your feedback.
Answers.UserDisplayName: snowflake
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T13:58:06.230Z
Answers.Text: Welcome to 'traditional' statistics @snowflake ... Q-Q is seen as a general approach. If it is approximately linear you good to go. In this instance that isn't quite the case. HW has meaning here (biological meaning). PCA doesn't have much meaning, its okay for data-mining, whilst GLM or GLMM is the sort of analysis where transformation is used to linearise a Q-Q. Peronally I'd stick with what you have and look at the biologial understanding. Just my opinion.
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T14:00:39.940Z
Answers.Text: @gringer post corrected
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T14:50:34.047Z
Answers.Text: @Michael Ok thanks, well that's a bit strange for me but I have to live with it :) But I appreciate your explanation and opinion a lot!
Answers.UserDisplayName: snowflake
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T16:05:37.790Z
Answers.Text: Hi @snowflake keep in mind the central rule of ML is that ANY equation can be replicated in ML, particularly deep learning. I'm not sure a population genecist will agree :-) but thats the theory. So you could use the classification to train an ML model and ditch the 'traditional approach' in future.
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T19:10:48.150Z
Answers.Text: PCA *are* used in GWAS; they're commonly used for covariate adjustment in the statistical models.
Answers.UserDisplayName: gringer
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-22T06:52:12.690Z
Answers.Text: @Michael That's the main goal and my interest. However I first need a running GWAS for comparison. I think looking at the entire SNP pattern instead of single SNPs  could be more informative (similar to  image patterns or finger prints). Let's see...
Answers.UserDisplayName: snowflake
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-22T12:05:14.603Z
Answers.Text: @snowflake okay I see, personally I think your HW results are fine given the caveats. However, if population genomic structure through HW isn't your thing (I would say it is a clear biological model) then first I'd do PCA and then depending on your purpose assess the populations e.g. via Fst. Personally I would just do a phylogenetic tree given these are SNPs. From a phylogenetics point of view the problem with PCA is that it is biologically difficult to interpret why the populations or different regions of the genome are separating. The thing which supercedes this dogma is ML partic. DNN
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-22T12:45:34.170Z
Answers.Text: A phylogenetic tree depends on linkage disequilibrium, which is the predominant basis of your data. So you could use your HW analysis to justify a phylogeny. Obviously I don't know the species.
Answers.UserDisplayName: M__
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-10-08T07:10:30.820Z
Answers.Text: Follow up question https://bioinformatics.stackexchange.com/questions/14500/gwas-phenotype-data-format-and-preprocessing
Answers.UserDisplayName: snowflake
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-20T18:07:55.287Z
Answers.LastActivityDate: 2020-09-21T13:59:19.157Z
Answers.LastEditDate: 2020-09-21T13:59:19.157Z
Answers.UserDisplayName: M__
Answers.Body: <p>Based on the QQ plots, filter out variants with p &gt; 0.01, then display the beta statistic (or whatever other statistic you're using to demonstrate association). If you want to incorporate the p-value in the displayed results, show the minimum test statistic from a 95% confidence interval. Here's an attempt I made at doing that with some GWAS summary statistics relating to household income (something that wouldn't be expected to have a heritable component):</p>
  <p><a href="https://twitter.com/gringene_bio/status/1207617723586371584" rel="nofollow noreferrer">https://twitter.com/gringene_bio/status/1207617723586371584</a></p>
  <p>After looking at the test statistic, you may be able to make better decisions about thresholding with the other values. One thing you could try is to incorporate those values into a scatter plot that shows more than two dimensions. The other statistics could be represented as size, colour, or Z location on a 3D scatter plot.</p>
  <p>P-values should not be overemphasised, used for ranking, or displayed on their own. See the American Statistical Association's report on the p-value for more information:</p>
  <p><a href="https://amstat.tandfonline.com/doi/full/10.1080/00031305.2016.1154108#_i30" rel="nofollow noreferrer">https://amstat.tandfonline.com/doi/full/10.1080/00031305.2016.1154108#_i30</a></p>
  <p>As far as I'm aware, the p-values as calculated by GWAS software represent uncertainty in the test statistic (i.e. t, or beta), and <em>not</em> the probability of association. If you want to calculate that association probability, carry out bootstrap sub-sampling of the individuals and then determine the proportion of sub-samples in which the tested variant appears in, say, the highest 5% of variants, ranked by the test statistic. Instead of filtering on HWE, MAF and MCR, you could try filtering on the bootstrap proportion instead, which should incorporate all statistics that introduce uncertainty into the dataset.</p>
  <p>More details and figures can be found on this poster:</p>
  <p><a href="https://f1000research.com/posters/5-2190" rel="nofollow noreferrer">https://f1000research.com/posters/5-2190</a></p>
  <p>If you're going to be generating a PCA on these data, please consider labelling the cases and controls (or low/high for quantitative traits) on the PCA to verify that the groups look similar in terms of overall population structure. Hopefully you won't find any blobs in the PCA that have only one of the two groups in them, but if you do, they should be dealt with first. Any differences in population structure (even very small, subtle differences) will dominate the GWAS results because structure-associated signal will cover a huge region of the genome in comparison to effects from a few (or a few hundred) SNPs. Here are some papers I've found insightful regarding this issue:</p>
  <ul>
  <li><a href="https://doi.org/10.1006/tpbi.2001.1543" rel="nofollow noreferrer">https://doi.org/10.1006/tpbi.2001.1543</a></li>
  <li><a href="https://doi.org/10.1038/s41467-018-05257-7" rel="nofollow noreferrer">https://doi.org/10.1038/s41467-018-05257-7</a></li>
  </ul>
  
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-21T06:49:39.633Z
Answers.Text: Thanks a lot for this detailed explanation, I will have a look on the beta statistics as well. plink is reporting beta and std, you are right.
Answers.UserDisplayName: snowflake
Answers.ContentLicense: CC BY-SA 4.0
Answers.CreationDate: 2020-09-19T20:01:04.687Z
Answers.LastActivityDate: 2020-09-21T11:08:19.177Z
Answers.LastEditDate: 2020-09-21T11:08:19.177Z
Answers.UserDisplayName: gringer
Body: <p>I did a bunch of GWAS analysis (linear model without covariates) with applying different quality controls. How to choose the optimal settings when filtering for minor allele frequency (maf), Hardy-Weinberg equilibrium (HW) and missing call rates (mcr)? Are there any rules? Here are some examples:</p>
  <p>Mahnhatten and QQ plot with settings maf=0.02, HW=1e-10, mcr = 0.2:</p>
  <p><a href="https://i.sstatic.net/YEgXO.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/YEgXO.png" alt="maf=0.02, HW=1e-10, missing call = 0.2" /></a></p>
  <p><a href="https://i.sstatic.net/lLKB7.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/lLKB7.png" alt="enter image description here" /></a></p>
  <p>Manhattan and QQ plot with settings maf=0.01, HW=0.0001, mcr = 0.01:</p>
  <p><a href="https://i.sstatic.net/Co5r7.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/Co5r7.png" alt="enter image description here" /></a></p>
  <p><a href="https://i.sstatic.net/shWXk.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/shWXk.png" alt="enter image description here" /></a></p>
  <p>Mahnhatten and QQ plot with settings maf=0.02, HW=1e-10, mcr = 0.05:
  <a href="https://i.sstatic.net/eqvnc.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/eqvnc.png" alt="enter image description here" /></a></p>
  <p><a href="https://i.sstatic.net/eqvnc.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/eqvnc.png" alt="enter image description here" /></a></p>
  <p><a href="https://i.sstatic.net/lOrHz.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/lOrHz.png" alt="enter image description here" /></a></p>
  
Comments.ContentLicense: CC BY-SA 4.0
Comments.CreationDate: 2020-09-19T15:59:57.223Z
Comments.Text: I don't think there is a principled way of choosing such cutoffs. The first and last set of filters look much better than the second set.
Comments.UserDisplayName: winni2k
Tags: <python><gwas><plink>
Title: Interpreting GWAS results with different settings
ThreadText: # Interpreting GWAS results with different settings
  
  (Asked by: username_0 on 2020-09-19)
  
  I did a bunch of GWAS analysis (linear model without covariates) with applying different quality controls. How to choose the optimal settings when filtering for minor allele frequency (maf), Hardy-Weinberg equilibrium (HW) and missing call rates (mcr)? Are there any rules? Here are some examples:
  
  Mahnhatten and QQ plot with settings maf=0.02, HW=1e-10, mcr = 0.2:
  
  [![maf=0.02, HW=1e-10, missing call = 0.2](https://i.sstatic.net/YEgXO.png)](https://i.sstatic.net/YEgXO.png)
  
  [![enter image description here](https://i.sstatic.net/lLKB7.png)](https://i.sstatic.net/lLKB7.png)
  
  Manhattan and QQ plot with settings maf=0.01, HW=0.0001, mcr = 0.01:
  
  [![enter image description here](https://i.sstatic.net/Co5r7.png)](https://i.sstatic.net/Co5r7.png)
  
  [![enter image description here](https://i.sstatic.net/shWXk.png)](https://i.sstatic.net/shWXk.png)
  
  Mahnhatten and QQ plot with settings maf=0.02, HW=1e-10, mcr = 0.05:
  [![enter image description here](https://i.sstatic.net/eqvnc.png)](https://i.sstatic.net/eqvnc.png)
  
  [![enter image description here](https://i.sstatic.net/eqvnc.png)](https://i.sstatic.net/eqvnc.png)
  
  [![enter image description here](https://i.sstatic.net/lOrHz.png)](https://i.sstatic.net/lOrHz.png)
  
  **Tags:** python, gwas, plink
  
  > **Comments:**
  > **username_1**: I don't think there is a principled way of choosing such cutoffs. The first and last set of filters look much better than the second set.  
  
  ---
  
  ## 2 answers
  
  ### Answer by username_2 on 2020-09-20. Score: 0 | Accepted answer
  
  The genomes under analysis are under predominant linkage disequilibrium, i.e. linkage. Chromosome 1 shows strong (or stronger - cut-off depending) HW signal as do one or two telomeres - telomeres are known to undergo recombination, so this may not be unsurprising. Nothing wrong with the HW analysis for the core conclusions, its the stuff between the red and blue lines thats difficult (Q-Q plots).
  
  The Q-Q plots are a bit more complex and I would go with the lower graph you presented together with its HW-analysis. There is clearly part of the data you are not accounting for, @username_3 says the Q-Q plot conform to the beta distribution .... they should be bang linear. Lets look in more detail,
  
  * 0 to 2 ... predominant disequilibrium that looks fine, nice and linear
  * 5.5+ .... HW signal that looks ok, nice and linear
  * 2 to 5.5 ... complicated
  
  Given you have a cut-off (red line) at around 5.5 then I'd sort of try and talk your way out of this one. You can do a full investigation, such as transformation to linearise the entire Q-Q plot. Personally I've never done this within HW analysis so I dunno, I've done it for GLM. I'd just say they between the blue and red line represents a complex HW signal which we can't be interpreted but either side of that boundary we're good to go. In other words keep it simple.
  
  > **Comments:**
  > **username_3**: "username_3 says the Q-Q plot conform to the beta distribution" - no, I didn't say this. The *p-values* should conform to the linear distribution if the association is not omnigenic, but the *test statistic* [which is not a p-value] will not fit a linear distribution. Most GWAS programs produce beta values as their association statistic, which is why I mentioned beta.  
  > **username_0**: Thanks Michael for your opinion. I am coming from the ML field were evaluation of a model with training data is more clear to me. I found it very hard to estimate GWAS filtered by "feeling" and looking at plots (also because I want to do it automatically). I am new to GWAS and all I need is a simple GWAS. I was unsure about using PCA or LMM  or I could go for a simple linear analysis and whether this looks kind of "normal".Thanks for your feedback.  
  > **username_2**: Welcome to 'traditional' statistics @username_0 ... Q-Q is seen as a general approach. If it is approximately linear you good to go. In this instance that isn't quite the case. HW has meaning here (biological meaning). PCA doesn't have much meaning, its okay for data-mining, whilst GLM or GLMM is the sort of analysis where transformation is used to linearise a Q-Q. Peronally I'd stick with what you have and look at the biologial understanding. Just my opinion.  
  > **username_2**: @username_3 post corrected  
  > **username_0**: @Michael Ok thanks, well that's a bit strange for me but I have to live with it :) But I appreciate your explanation and opinion a lot!  
  > **username_2**: Hi @username_0 keep in mind the central rule of ML is that ANY equation can be replicated in ML, particularly deep learning. I'm not sure a population genecist will agree :-) but thats the theory. So you could use the classification to train an ML model and ditch the 'traditional approach' in future.  
  > **username_3**: PCA *are* used in GWAS; they're commonly used for covariate adjustment in the statistical models.  
  > **username_0**: @Michael That's the main goal and my interest. However I first need a running GWAS for comparison. I think looking at the entire SNP pattern instead of single SNPs  could be more informative (similar to  image patterns or finger prints). Let's see...  
  > **username_2**: @username_0 okay I see, personally I think your HW results are fine given the caveats. However, if population genomic structure through HW isn't your thing (I would say it is a clear biological model) then first I'd do PCA and then depending on your purpose assess the populations e.g. via Fst. Personally I would just do a phylogenetic tree given these are SNPs. From a phylogenetics point of view the problem with PCA is that it is biologically difficult to interpret why the populations or different regions of the genome are separating. The thing which supercedes this dogma is ML partic. DNN  
  > **username_2**: A phylogenetic tree depends on linkage disequilibrium, which is the predominant basis of your data. So you could use your HW analysis to justify a phylogeny. Obviously I don't know the species.  
  > **username_0**: Follow up question https://bioinformatics.stackexchange.com/questions/14500/gwas-phenotype-data-format-and-preprocessing  
  
  ---
  
  ### Answer by username_3 on 2020-09-19. Score: 2
  
  Based on the QQ plots, filter out variants with p > 0.01, then display the beta statistic (or whatever other statistic you're using to demonstrate association). If you want to incorporate the p-value in the displayed results, show the minimum test statistic from a 95% confidence interval. Here's an attempt I made at doing that with some GWAS summary statistics relating to household income (something that wouldn't be expected to have a heritable component):
  
  <https://twitter.com/gringene_bio/status/1207617723586371584>
  
  After looking at the test statistic, you may be able to make better decisions about thresholding with the other values. One thing you could try is to incorporate those values into a scatter plot that shows more than two dimensions. The other statistics could be represented as size, colour, or Z location on a 3D scatter plot.
  
  P-values should not be overemphasised, used for ranking, or displayed on their own. See the American Statistical Association's report on the p-value for more information:
  
  <https://amstat.tandfonline.com/doi/full/10.1080/00031305.2016.1154108#_i30>
  
  As far as I'm aware, the p-values as calculated by GWAS software represent uncertainty in the test statistic (i.e. t, or beta), and *not* the probability of association. If you want to calculate that association probability, carry out bootstrap sub-sampling of the individuals and then determine the proportion of sub-samples in which the tested variant appears in, say, the highest 5% of variants, ranked by the test statistic. Instead of filtering on HWE, MAF and MCR, you could try filtering on the bootstrap proportion instead, which should incorporate all statistics that introduce uncertainty into the dataset.
  
  More details and figures can be found on this poster:
  
  <https://f1000research.com/posters/5-2190>
  
  If you're going to be generating a PCA on these data, please consider labelling the cases and controls (or low/high for quantitative traits) on the PCA to verify that the groups look similar in terms of overall population structure. Hopefully you won't find any blobs in the PCA that have only one of the two groups in them, but if you do, they should be dealt with first. Any differences in population structure (even very small, subtle differences) will dominate the GWAS results because structure-associated signal will cover a huge region of the genome in comparison to effects from a few (or a few hundred) SNPs. Here are some papers I've found insightful regarding this issue:
  
  * <https://doi.org/10.1006/tpbi.2001.1543>
  * <https://doi.org/10.1038/s41467-018-05257-7>
  
  > **Comments:**
  > **username_0**: Thanks a lot for this detailed explanation, I will have a look on the beta statistics as well. plink is reporting beta and std, you are right.  
  
  ---
  
  