difference in mapping genes w/ different MHC analyses

[Isotta Bozzarelli]

Hello,

I have a doubt about the results I obtain from the analyses of one of my GWAS that have a hit in the MHC region. I tried to run it 

1) with the MHC region; 

2)w/out the MHC region only for MAGMA; 

3) w/out the MHC region only for annotation; 

4) w/out for both.

for the first case i obtain 461 mapped genes, while I get 359 mapped genes in the second case. this is weird to me because I only excluded the MHC region from MAGMA and not from annotation. I checked the ensg of the missing genes and they all belong to chr6. 

in the 3rd and 4th case i only have 211 genes, but checking the missing ensg I found out that these are genes missing not only from chr6 but from all chromosomes. 

Could you explain to me why I see these differences?

thanks, regards,

Isotta

[Tanya Phung]

Hi Isotta, 

Could you please list all of the jobIDs? 

Best,
Tanya

[Isotta Bozzarelli]

Hello Tanya,
the jobs are:
1) 588876 (with MHC)
2) 588882 (no MHC in magma)
3) 588881  (no MHC in annotation)

4) 588877 (without MHC)

I have also tried to do other attempts to understand better how it works but I get more and more confused. I performed the same analysis including MHC for all but not performing MAGMA (analysis number 589536, so not ticking MAGMA analysis box) and I get the same number of mapped genes as 588876 (with MHC, and running magma, 461 genes). Then I performed the analysis without MHC ONLY FOR MAGMA, and NOT RUNNING magma (not ticking it) (analysis number 589551) and I have the same results as in 588882 (no MHC in magma), so 359 mapped genes. This is very confusing to me.

Does MAGMA analysis affects in any way the mapping of the genes? From the tutorial I have understood that it shouldn't (point 3, so mapping, is before and separated from MAGMA) is my understanding correct?
If it doesn't affect it, then why I get different results in my mapping genes when I select MHC included vs MHC excluded only for MAGMA?
If it affects it, then why I get the same results also when I select the same two parameters as above for MHC but I select of not performing MAGMA analysis?
Thanks for your time,
Best,
Isotta

[Tanya Phung]

Hi Isotta, 

First, thank you for providing the details. 

Second, this issue appears to be a bug in the code base. I will try to explain below.

• You are correct that the expected behavior should be: 
  • When you select to exclude MHC in MAGMA analyses only, then the mapped genes from the annotation steps should not be affected. However, this is currently not the case because of a bug in the code. 
  • In job 588876, MHC is included and there were 461 genes. 
  • In job 588882, even though you exclude MHC in MAGMA only, currently what the code does is that: 
    • In the step to generate the genomic risk loci, it does this correctly, meaning MHC is not excluded. You can verify this by the fact that they have the same genomic risk loci
      • cat 588876/GenomicRiskLoci.txt | wc -l
        21 #20 genomic risk loci because the first line is the header 
        
        cat 588882/GenomicRiskLoci.txt | wc -l
        21  #20 genomic risk loci because the first line is the header 
        
        diff 588876/GenomicRiskLoci.txt 588882/GenomicRiskLoci.txt (returns nothing because they are the same)
    • However, at the gene mapping step, restricting MHC exclusion to MAGMA only is not taken into account and the code in the gene mapping step incorrectly removes MHC. This results in no genes falling in the MHC region being mapped. 
  • In job 588881 (and similarly in job 588877), when MHC is excluded in the annotation step, MHC is excluded in the step to generate the genomic risk loci. Therefore, in the previous 2 jobs (588876 and  588882) there were 20 genomic loci but for jobs  588881 and 588877 there are 19 genomic risk loci because the genomic risk loci that fall within the MHC region is excluded. You can verify this yourself. 
    • Because 1 genomic risk loci (6:31327400:C:T) is removed. Subsequently in the gene mapping step, no genes are mapped to this genomic risk loci. 
    • Now, why is the number of mapped genes 211 genes (while there are 359 mapped genes in job 588882). This is because of the eQTL mapping. 
      • Genomic risk locus #9 (6:31327400:C:T  rs2523587       6       31327400) from job 588882 has the following mapped genes by eQLT mapping. Below is a few of the genes that are mapped to this genomic loci by eQTL mapping. The last column (red squares) indicate that this is genomic risk loci #9. Because this is eQTL mapping and not positional mapping, the mapped genes are not necessarily on chr6. 
      • In job 588882 these eQTLs mapped genes are not excluded because they are mapped outside of the MHC region. 

I hope that my explanation makes sense. In short, the reason why you are observing these discrepancies is because in the gene mapping step, there is a bug in the code where it is not recognizing that excluding MHC is restricting to MAGMA only. 

Third, because this is a bug. I will need to fix the code. Please expect 3-4 weeks (maybe more) for a bug fix. Thank you for your patience. 

Best,
Tanya

[Tanya Phung]

Hi Isotta, 

This issue is now fixed: https://groups.google.com/g/fuma-gwas-users/c/dSqfBYXz9Wg

I checked this using the example input data (and GTEx v8 Heart data for eQTL mapping): 

1. MHC is included:

   1. https://fuma.ctglab.nl/browse/603621#summaryTable

   2. 326 mapped genes

2. MHC is excluded in MAGMA only

   2. 326 mapped genes

3. MHC is excluded in annotation only

   1. https://fuma.ctglab.nl/browse/603623#summaryTable

   2. 278 mapped genes

Best,
Tanya