How to determine a significant WAIC threshold in a GWAS setting?

[Xulong Wang]

Hi Stan,

I just implemented my model with Stan in a genome-wide association study. It worked great. Detailedly, I have 1 million SNP. I run my model for each of them independently. With that, I computed WAIC to quantify the fitness. Now I have 1 million WAIC quantities. I also have a null model, in which no SNP was included. By subtracting the WAIC quantity of the null model, I have the WAIC differences for each of the SNP, which describes how much the model was improved by including each of the SNP.

How to determine a cutoff for the WAIC differences (or original WAIC) to identify the true association? This is reminiscent of multiple comparison correction of p-value in the frequentist framework.

Can I generate a null distribution of the WAIC by shuffling the genotype of each SNP (permutation), from which to compute a p-value for each of the real WAIC?

Also, given that the full posterior distribution of each SNP's effect size was obtained by sampling, can I simply use the 95% CI and not worry about the multiple comparison issue at all?

Xulong

[Andrew Gelman]

Hi, if you are comparing just 2 models, you can get a s.e. for Waic as described in this paper:

http://www.stat.columbia.edu/~gelman/research/unpublished/loo_stan.pdf

With more than 2 models to compare, it would make sense to fit some sort of hierarchical model to the Waic output; we’ve discussed this idea but have never done it.

A

--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Aki Vehtari]

On Monday, November 30, 2015 at 10:23:29 PM UTC+2, Xulong Wang wrote:

> Can I generate a null distribution of the WAIC by shuffling the genotype of each SNP (permutation), from which to compute a p-value for each of the real WAIC?

Yes, you can. This way will get some idea how big difference in WAIC is likely to be significant, but making independent models for each SNP is bad idea.

> Also, given that the full posterior distribution of each SNP's effect size was obtained by sampling, can I simply use the 95% CI and not worry about the multiple comparison issue at all?

If you make independent models, you need to worry about the multiple comparison. 

It would be better to make one hierarchical model for all SNPs. Hierarchcial model with 1 million SNPs is challenging, but here are couple papers

http://dx.doi.org/10.1371/journal.pone.0049445

http://dx.doi.org/10.1371/journal.pone.0029115

and code

http://becs.aalto.fi/en/research/bayes/bmagwa/

(unfortunately this is not possible yet with Stan)

If the hierarchical model is good, then you don't need to worry about the multiple comparison (see on Andrew's papers).

On Tuesday, December 1, 2015 at 7:17:13 AM UTC+2, Andrew Gelman wrote:

> With more than 2 models to compare, it would make sense to fit some sort of hierarchical model to the Waic output; we’ve discussed this idea but have never done it.

It's better to make hierarchical model for SNPs than for WAIC output. Making a hierarchical model for WAIC output would be as difficult as for SNP data, but less efficient as you need to make first the independent models, and WAIC summary will discard information compared to the original data.

Aki

[Andrew Gelman]

> On Dec 1, 2015, at 6:55 AM, Aki Vehtari <Aki.V...@aalto.fi> wrote:
>
> It's better to make hierarchical model for SNPs than for WAIC output. Making a hierarchical model for WAIC output would be as difficult as for SNP data, but less efficient as you need to make first the independent models, and WAIC summary will discard information compared to the original data.

Yes, I agree.

[Xulong Wang]

Andrew and Aki, 

Thanks much for the precious advices.

Hierarchical modeling for GWAS in Bayesian is very appealing to me. I will try to implement this method in my study. I will come back and report my progress.

Best,

Xulong

[Aki Vehtari]

On Wednesday, December 2, 2015 at 6:23:18 PM UTC+2, Xulong Wang wrote:

Hierarchical modeling for GWAS in Bayesian is very appealing to me. I will try to implement this method in my study. I will come back and report my progress.

Which method you will try to implement? There are already implementations of hierarchical modeling for GWAS in Bayesian way, e.g. this one

http://becs.aalto.fi/en/research/bayes/bmagwa/

so you could start by testing those, before implementing something new.

Aki

[Xulong Wang]

Aki,

Thanks for the advice. I will certainly start from bmagwa, to get familiar with the use of hierarchical modeling for GWA. But, if not too overwhelming, I might eventually need to build something new, because my model is an ordered categorical mixed model with known covariate matrix.

[Xulong Wang]

Andrew, with s.e of WAIC statistics available, would it be reasonable to estimate a p-value by a normal distribution? 

[Andrew Gelman]

Yes, I think so.

A