JonJala / mtag

Python command line tool for Multi-Trait Analysis of GWAS (MTAG)
GNU General Public License v3.0
169 stars 54 forks source link

High FDR (41%) when MTAGing 4 traits but pairwise MTAG always <5% #147

Open pjordab opened 2 years ago

pjordab commented 2 years ago

Hi users and developers,

Many thanks in advance for your time and help.

I am computing MTAG for four traits and, surprisingly, I get an FDR of 41%. However, the FDRs of each of the traits in the pairwise anaysis is < 5%.

One possible explanation is the high chi-square of trait 4 and its low correlation with the other traits, but I do not understand why this only affects the analysis with 4 traits and not the paired analysis. I do not know if you know what could be the explanation.

See here a summary of the MTAG results, the pairwise FDRs and rg, and the mean chi 2 for each trait.

SNPs FDR Chi 2 TRAIT 1 pre-MTAG Chi 2 TRAIT 1 - post MTAG Maximum N 5883227 0.41 1.138 1.317 2279537

Pairwise FDR (always < 0.05)

     TRAIT 1    TRAIT 2 TRAIT 3 TRAIT 4

TRAIT 1 - 0.020 0.028 0.040 TRAIT 2 0.004 - 0.007 0.005 TRAIT 3 0.003 0.003 - 0.002 TRAIT 4 7.10E-06 7.20E-06 7.39E-06 -

RG (low for TRAIT 4, and for TRAIT 2 ~TRAIT 3)

           TRAIT 1  TRAIT 2 TRAIT 3 TRAIT 4

TRAIT 1 - 0.71 0.56 0.28 TRAIT 2 0.71 - 0.18 0.33 TRAIT 3 0.56 0.18 - 0.12 TRAIT 4 0.28 0.33 0.12 -

Chi squared (TRAIT 4 sumstats have higher chi squared compared to the other sumstats)

TRAIT 1 1.138 TRAIT 2 1.078 TRAIT 3 1.03 TRAIT 4 2.838

Many many thanks for your help!

paturley commented 2 years ago

Hello,

I'm surprised that you even got the maxFDR code to run for 4 traits. It can get quite slow even with 3. What options did you use for the 4 trait version vs the 2 trait version?

Best, Patrick

On Fri, Nov 26, 2021 at 11:20 AM pjordab @.***> wrote:

Hi users and developers,

Many thanks in advance for your time and help.

I am computing MTAG for four traits and, surprisingly, I get an FDR of 41%. However, the FDRs of each of the traits in the pairwise anaysis is < 5%.

One possible explanation is the high chi-square of trait 4 and its low correlation with the other traits, but I do not understand why this only affects the analysis with 4 traits and not the paired analysis. I do not know if you know what could be the explanation.

See here a summary of the MTAG results, the pairwise FDRs and rg, and the mean chi 2 for each trait.

SNPs FDR Chi 2 TRAIT 1 pre-MTAG Chi 2 TRAIT 1 - post MTAG Maximum N 5883227 0.41 1.138 1.317 2279537

Pairwise FDR (always < 0.05)

                              Second trait
             TRAIT 1  TRAIT 2 TRAIT 3 TRAIT 4

First trait TRAIT 1 - 0.020 0.028 0.040 TRAIT 2 0.004 - 0.007 0.005 TRAIT 3 0.003 0.003 - 0.002 TRAIT 4 7.10E-06 7.20E-06 7.39E-06 -

RG (low for TRAIT 4, and for TRAIT 2 ~TRAIT 3)

                                    Second trait
             TRAIT 1  TRAIT 2 TRAIT 3 TRAIT 4

First trait TRAIT 1 - 0.71 0.56 0.28 TRAIT 2 0.71 - 0.18 0.33 TRAIT 3 0.56 0.18 - 0.12 TRAIT 4 0.28 0.33 0.12 -

Chi squared (TRAIT 4 sumstats have higher chi squared compared to the other sumstats)

Chi 2 TRAIT 1 1.138 TRAIT 2 1.078 TRAIT 3 1.03 TRAIT 4 2.838

Many many thanks for your help!

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/JonJala/mtag/issues/147, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFBUB5P4AK4SPR6JHHZWT7DUN6XTRANCNFSM5I27SLCQ . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

pjordab commented 2 years ago

Hi,

The same ones for the 2 traits and 4 traits analysis:

Calling ./mtag.py \ --stream-stdout \ --n-min 0.0 \ --sumstats trait1,trait4 \ --fdr \ --out ./MTAG_nov2021_pairwise-fdr/trait1-trait4

Thank you! Paloma

paturley commented 2 years ago

Do you have the log file for the 4-trait maxFDR run and one of the 2-trait maxFDR runs?

On Wed, Dec 1, 2021 at 9:43 AM pjordab @.***> wrote:

Hi,

The same ones for the 2 traits and 4 traits analysis:

Calling ./mtag.py --stream-stdout --n-min 0.0 --sumstats trait1,trait4 --fdr --out ./MTAG_nov2021_pairwise-fdr/trait1-trait4

Thank you! Paloma

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/JonJala/mtag/issues/147#issuecomment-983709611, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFBUB5JPQENRVCTRPG3UFUTUOYYCFANCNFSM5I27SLCQ . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

pjordab commented 2 years ago

Hi Patrick and thanks for your help, Researching of possible for this FDR I found a post in which it indicates that for binary phenotypes (from case-control studies) I need to calculate the effective n instead of the sample size.

https://github.com/JonJala/mtag/issues/15

https://github.com/JonJala/mtag/issues/10

I'll do it and come back to you with my new results.

Many thanks!

Paloma

pjordab commented 2 years ago

LOG FILE FDR 41% - 4 TRAIT ANALYSIS

2021/12/16/03:47:14 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> <> <> MTAG: Multi-trait Analysis of GWAS <> Version: 1.0.8 <> (C) 2017 Omeed Maghzian, Raymond Walters, and Patrick Turley <> Harvard University Department of Economics / Broad Institute of MIT and Harvard <> GNU General Public License v3 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> <> Note: It is recommended to run your own QC on the input before using this program. <> Software-related correspondence: maghzian@nber.org <> All other correspondence: paturley@broadinstitute.org <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

Calling ./mtag.py \ --stream-stdout \ --sumstats TRAIT_1,TRAIT_2,TRAIT_3,TRAIT_4 \ --fdr \ --out ./ANALYSIS

2021/12/16/03:47:14 PM Beginning MTAG analysis... 2021/12/16/03:47:14 PM MTAG will use the Z column for analyses. 2021/12/16/03:47:25 PM Read in Trait 1 summary statistics (6967088 SNPs) from TRAIT_1 ... 2021/12/16/03:47:25 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:47:25 PM Munging Trait 1 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/03:47:25 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:47:25 PM Interpreting column names as follows: 2021/12/16/03:47:25 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/03:47:26 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/03:47:36 PM Read 6967088 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 0 SNPs with out-of-bounds p-values. Removed 0 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 6967088 SNPs remain. 2021/12/16/03:47:41 PM Removed 0 SNPs with duplicated rs numbers (6967088 SNPs remain). 2021/12/16/03:47:42 PM Removed 0 SNPs with N < 120050.470727 (6967088 SNPs remain). 2021/12/16/03:48:53 PM Median value of SIGNED_SUMSTAT was 0.00666666666667, which seems sensible. 2021/12/16/03:48:54 PM Dropping snps with null values 2021/12/16/03:48:54 PM Metadata: 2021/12/16/03:48:55 PM Mean chi^2 = 1.143 2021/12/16/03:48:55 PM Lambda GC = 1.114 2021/12/16/03:48:55 PM Max chi^2 = 83.718 2021/12/16/03:48:55 PM 298 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/03:48:55 PM Conversion finished at Thu Dec 16 15:48:55 2021 2021/12/16/03:48:55 PM Total time elapsed: 1.0m:29.54s 2021/12/16/03:49:06 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:49:06 PM Munging of Trait 1 complete. SNPs remaining: 6967088 2021/12/16/03:49:06 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

2021/12/16/03:49:33 PM Read in Trait 2 summary statistics (8557337 SNPs) from TRAIT_2 ... 2021/12/16/03:49:33 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:49:33 PM Munging Trait 2 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/03:49:33 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:49:33 PM Interpreting column names as follows: 2021/12/16/03:49:33 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/03:49:33 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/03:49:45 PM Read 8557337 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 0 SNPs with out-of-bounds p-values. Removed 551192 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 8006145 SNPs remain. 2021/12/16/03:49:51 PM Removed 0 SNPs with duplicated rs numbers (8006145 SNPs remain). 2021/12/16/03:49:52 PM Removed 0 SNPs with N < 149818.257102 (8006145 SNPs remain). 2021/12/16/03:51:12 PM Median value of SIGNED_SUMSTAT was -0.00122699386503, which seems sensible. 2021/12/16/03:51:12 PM Dropping snps with null values 2021/12/16/03:51:13 PM Metadata: 2021/12/16/03:51:13 PM Mean chi^2 = 1.088 2021/12/16/03:51:14 PM Lambda GC = 1.007 2021/12/16/03:51:14 PM Max chi^2 = 458.046 2021/12/16/03:51:14 PM 2610 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/03:51:14 PM Conversion finished at Thu Dec 16 15:51:14 2021 2021/12/16/03:51:14 PM Total time elapsed: 1.0m:40.94s 2021/12/16/03:51:27 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:51:27 PM Munging of Trait 2 complete. SNPs remaining: 8006145 2021/12/16/03:51:27 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

2021/12/16/03:51:57 PM Read in Trait 3 summary statistics (9397905 SNPs) from TRAIT_3 ... 2021/12/16/03:51:57 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:51:57 PM Munging Trait 3 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/03:51:57 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:51:57 PM Interpreting column names as follows: 2021/12/16/03:51:57 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/03:51:57 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/03:52:02 PM WARNING: 111 SNPs had P outside of (0,1]. The P column may be mislabeled. 2021/12/16/03:52:11 PM Read 9397905 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 111 SNPs with out-of-bounds p-values. Removed 607579 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 8790215 SNPs remain. 2021/12/16/03:52:16 PM Removed 5619 SNPs with duplicated rs numbers (8784596 SNPs remain). 2021/12/16/03:52:17 PM Removed 0 SNPs with N < 152147.044198 (8784596 SNPs remain). 2021/12/16/03:53:48 PM Median value of SIGNED_SUMSTAT was -0.00709219858156, which seems sensible. 2021/12/16/03:53:49 PM Dropping snps with null values 2021/12/16/03:53:49 PM Metadata: 2021/12/16/03:53:50 PM Mean chi^2 = 1.376 2021/12/16/03:53:50 PM Lambda GC = 1.215 2021/12/16/03:53:50 PM Max chi^2 = 1452.322 2021/12/16/03:53:50 PM 12398 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/03:53:50 PM Conversion finished at Thu Dec 16 15:53:50 2021 2021/12/16/03:53:50 PM Total time elapsed: 1.0m:53.48s 2021/12/16/03:54:04 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:54:04 PM Munging of Trait 3 complete. SNPs remaining: 8790215 2021/12/16/03:54:04 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

2021/12/16/03:54:20 PM Trait 3: Dropped 5619 SNPs for duplicate values in the "snp_name" column 2021/12/16/03:54:30 PM Read in Trait 4 summary statistics (7087911 SNPs) from TRAIT_4 ... 2021/12/16/03:54:30 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:54:30 PM Munging Trait 4 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/03:54:30 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:54:30 PM Interpreting column names as follows: 2021/12/16/03:54:30 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/03:54:31 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/03:54:41 PM Read 7087911 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 0 SNPs with out-of-bounds p-values. Removed 0 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 7087911 SNPs remain. 2021/12/16/03:54:46 PM Removed 30 SNPs with duplicated rs numbers (7087881 SNPs remain). 2021/12/16/03:54:47 PM Removed 0 SNPs with N < 497212.0 (7087881 SNPs remain). 2021/12/16/03:55:58 PM Median value of SIGNED_SUMSTAT was 0.00967742, which seems sensible. 2021/12/16/03:55:58 PM Dropping snps with null values 2021/12/16/03:55:58 PM Metadata: 2021/12/16/03:55:59 PM Mean chi^2 = 2.838 2021/12/16/03:55:59 PM Lambda GC = 1.916 2021/12/16/03:55:59 PM Max chi^2 = 626.746 2021/12/16/03:55:59 PM 77637 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/03:55:59 PM Conversion finished at Thu Dec 16 15:55:59 2021 2021/12/16/03:55:59 PM Total time elapsed: 1.0m:28.72s 2021/12/16/03:56:11 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:56:11 PM Munging of Trait 4 complete. SNPs remaining: 7087911 2021/12/16/03:56:11 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/03:56:24 PM Trait 4: Dropped 30 SNPs for duplicate values in the "snp_name" column 2021/12/16/03:56:29 PM Dropped 1039218 SNPs due to strand ambiguity, 5927870 SNPs remain in intersection after merging trait1 2021/12/16/03:56:51 PM Flipped the signs of of 2904828 SNPs to make them consistent with the effect allele orderings of the first trait. 2021/12/16/03:56:55 PM Dropped 0 SNPs due to strand ambiguity, 5812916 SNPs remain in intersection after merging trait2 2021/12/16/03:57:18 PM Dropped 11 SNPs due to inconsistent allele pairs from phenotype 3. 5810375 SNPs remain. 2021/12/16/03:57:26 PM Flipped the signs of of 5810375 SNPs to make them consistent with the effect allele orderings of the first trait. 2021/12/16/03:57:31 PM Dropped 0 SNPs due to strand ambiguity, 5810375 SNPs remain in intersection after merging trait3 2021/12/16/03:57:59 PM Flipped the signs of of 3068650 SNPs to make them consistent with the effect allele orderings of the first trait. 2021/12/16/03:58:04 PM Dropped 0 SNPs due to strand ambiguity, 5617924 SNPs remain in intersection after merging trait4 2021/12/16/03:58:04 PM ... Merge of GWAS summary statistics complete. Number of SNPs: 5617924 2021/12/16/03:58:31 PM Using 5617924 SNPs to estimate Omega (0 SNPs excluded due to strand ambiguity) 2021/12/16/03:58:31 PM Estimating sigma.. 2021/12/16/04:01:32 PM Checking for positive definiteness .. 2021/12/16/04:01:32 PM Sigma hat: [[1.006 0.056 0.246 0.018] [0.056 0.853 0.051 0.028] [0.246 0.051 1.054 0.015] [0.018 0.028 0.015 1.075]] 2021/12/16/04:01:32 PM Beginning estimation of Omega ... 2021/12/16/04:01:33 PM Using GMM estimator of Omega .. 2021/12/16/04:01:35 PM Checking for positive definiteness .. 2021/12/16/04:01:35 PM Completed estimation of Omega ... 2021/12/16/04:01:35 PM Beginning MTAG calculations... 2021/12/16/04:01:55 PM ... Completed MTAG calculations. 2021/12/16/04:01:55 PM Writing Phenotype 1 to file ... 2021/12/16/04:02:48 PM Writing Phenotype 2 to file ... 2021/12/16/04:03:46 PM Writing Phenotype 3 to file ... 2021/12/16/04:04:43 PM Writing Phenotype 4 to file ... 2021/12/16/04:05:37 PM Summary of MTAG results:

Trait # SNPs used ... MTAG mean chi^2 GWAS equiv. (max) N 1 ...TRAIT_1 5617924 ... 1.323 423601 2 ...TRAIT_2 5617924 ... 1.412 308532 3 ...TRAIT_3 5617924 ... 1.382 232789 4 ...TRAIT_4 5617924 ... 2.699 754672

[4 rows x 7 columns]

Estimated Omega: [[7.682e-07 6.180e-07 6.299e-07 3.552e-07] [6.180e-07 1.143e-06 2.200e-07 5.683e-07] [6.299e-07 2.200e-07 1.731e-06 1.743e-07] [3.552e-07 5.683e-07 1.743e-07 2.452e-06]]

(Correlation): [[1. 0.66 0.546 0.259] [0.66 1. 0.156 0.339] [0.546 0.156 1. 0.085] [0.259 0.339 0.085 1. ]]

Estimated Sigma: [[1.006 0.056 0.246 0.018] [0.056 0.853 0.051 0.028] [0.246 0.051 1.054 0.015] [0.018 0.028 0.015 1.075]]

(Correlation): [[1. 0.06 0.238 0.017] [0.06 1. 0.053 0.029] [0.238 0.053 1. 0.014] [0.017 0.029 0.014 1. ]]

MTAG weight factors: (average across SNPs) [1.264 1.275 1.172 1.058]

2021/12/16/04:05:37 PM 2021/12/16/04:05:37 PM MTAG results saved to file. 2021/12/16/04:05:37 PM Beginning maxFDR calculations. Depending on the number of grid points specified, this might take some time... 2021/12/16/04:05:37 PM T=4 2021/12/16/04:22:21 PM Number of gridpoints to search: 80670 2021/12/16/04:22:21 PM Performing grid search using 1 cores. 2021/12/16/04:24:36 PM Grid search: 10.0 percent finished for . Time: 2.233 min 2021/12/16/04:26:52 PM Grid search: 20.0 percent finished for . Time: 4.508 min 2021/12/16/04:29:06 PM Grid search: 30.0 percent finished for . Time: 6.740 min 2021/12/16/04:31:22 PM Grid search: 40.0 percent finished for . Time: 9.000 min 2021/12/16/04:33:35 PM Grid search: 50.0 percent finished for . Time: 11.224 min 2021/12/16/04:35:48 PM Grid search: 60.0 percent finished for . Time: 13.444 min 2021/12/16/04:38:01 PM Grid search: 70.0 percent finished for . Time: 15.657 min 2021/12/16/04:40:16 PM Grid search: 80.0 percent finished for . Time: 17.896 min 2021/12/16/04:42:29 PM Grid search: 90.0 percent finished for . Time: 20.125 min 2021/12/16/04:44:43 PM Grid search: 100.0 percent finished for . Time: 22.353 min 2021/12/16/04:44:46 PM Saved calculations of fdr over grid points in ANALYSIS_fdr_mat.txt 2021/12/16/04:44:46 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/04:44:46 PM grid point indices for max FDR for each trait: [ 3782 5972 21444 42739] 2021/12/16/04:44:46 PM Maximum FDR 2021/12/16/04:44:46 PM Max FDR of Trait 1: 0.409029765019 at probs = [0. 0. 0. 0. 0. 0. 0. 0.1 0. 0. 0.4 0. 0.2 0.3 0. 0. ] 2021/12/16/04:44:46 PM Max FDR of Trait 2: 0.00929616247232 at probs = [0. 0. 0. 0. 0. 0.1 0. 0. 0. 0.1 0. 0. 0.4 0. 0.3 0.1] 2021/12/16/04:44:46 PM Max FDR of Trait 3: 0.0024796685601 at probs = [0. 0. 0.1 0.1 0. 0. 0. 0. 0. 0. 0.4 0. 0. 0.2 0.2 0. ] 2021/12/16/04:44:46 PM Max FDR of Trait 4: 6.94111066721e-06 at probs = [0. 0.3 0. 0. 0. 0.2 0. 0. 0. 0. 0. 0.1 0. 0.1 0.3 0. ] 2021/12/16/04:44:46 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/04:44:46 PM Completed FDR calculations. 2021/12/16/04:44:46 PM MTAG complete. Time elapsed: 57.0m:31.7294001579s

pjordab commented 2 years ago

LOG FILE - 2 TRAIT ANALYSIS

Please find attached the record of the pairwise analysis of TRAIT_1 and TRAIT_4 from previous analyses. TRAIT_4 is the one that has the highest mean chi 2 with respect to the other traits and TRAIT_1 and so I initially suspected that it might be the one that was causing an increase in the FDR in the TRAIT_1 mtag output.

2021/12/16/05:57:05 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> <> <> MTAG: Multi-trait Analysis of GWAS <> Version: 1.0.8 <> (C) 2017 Omeed Maghzian, Raymond Walters, and Patrick Turley <> Harvard University Department of Economics / Broad Institute of MIT and Harvard <> GNU General Public License v3 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> <> Note: It is recommended to run your own QC on the input before using this program. <> Software-related correspondence: maghzian@nber.org <> All other correspondence: paturley@broadinstitute.org <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

Calling ./mtag.py \ --stream-stdout \ --sumstats TRAIT_1,TRAIT_4 \ --fdr \ --out ./PAIRWISE_ANALYSIS

2021/12/16/05:57:05 PM Beginning MTAG analysis... 2021/12/16/05:57:05 PM MTAG will use the Z column for analyses. 2021/12/16/05:57:17 PM Read in Trait 1 summary statistics (6967088 SNPs) from TRAIT_1 ... 2021/12/16/05:57:17 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/05:57:17 PM Munging Trait 1 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/05:57:17 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/05:57:17 PM Interpreting column names as follows: 2021/12/16/05:57:17 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/05:57:17 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/05:57:27 PM Read 6967088 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 0 SNPs with out-of-bounds p-values. Removed 0 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 6967088 SNPs remain. 2021/12/16/05:57:32 PM Removed 0 SNPs with duplicated rs numbers (6967088 SNPs remain). 2021/12/16/05:57:33 PM Removed 0 SNPs with N < 120050.470727 (6967088 SNPs remain). 2021/12/16/05:58:46 PM Median value of SIGNED_SUMSTAT was 0.00666666666667, which seems sensible. 2021/12/16/05:58:46 PM Dropping snps with null values 2021/12/16/05:58:47 PM Metadata: 2021/12/16/05:58:47 PM Mean chi^2 = 1.143 2021/12/16/05:58:47 PM Lambda GC = 1.114 2021/12/16/05:58:48 PM Max chi^2 = 83.718 2021/12/16/05:58:48 PM 298 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/05:58:48 PM Conversion finished at Thu Dec 16 17:58:48 2021 2021/12/16/05:58:48 PM Total time elapsed: 1.0m:31.08s 2021/12/16/05:58:58 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/05:58:58 PM Munging of Trait 1 complete. SNPs remaining: 6967088 2021/12/16/05:58:58 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

2021/12/16/05:59:21 PM Read in Trait 2 summary statistics (7087911 SNPs) from TRAIT_4 ... 2021/12/16/05:59:21 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/05:59:21 PM Munging Trait 2 <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><>< 2021/12/16/05:59:21 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/05:59:21 PM Interpreting column names as follows: 2021/12/16/05:59:21 PM snpid: Variant ID (e.g., rs number) n: Sample size a1: a1, interpreted as ref allele for signed sumstat. pval: p-Value a2: a2, interpreted as non-ref allele for signed sumstat. z: Directional summary statistic as specified by --signed-sumstats.

2021/12/16/05:59:21 PM Reading sumstats from provided DataFrame into memory 10000000 SNPs at a time. 2021/12/16/05:59:31 PM Read 7087911 SNPs from --sumstats file. Removed 0 SNPs with missing values. Removed 0 SNPs with INFO <= None. Removed 0 SNPs with MAF <= 0.01. Removed 0 SNPs with SE <0 or NaN values. Removed 0 SNPs with out-of-bounds p-values. Removed 0 variants that were not SNPs. Note: strand ambiguous SNPs were not dropped. 7087911 SNPs remain. 2021/12/16/05:59:36 PM Removed 30 SNPs with duplicated rs numbers (7087881 SNPs remain). 2021/12/16/05:59:37 PM Removed 0 SNPs with N < 497212.0 (7087881 SNPs remain). 2021/12/16/06:00:48 PM Median value of SIGNED_SUMSTAT was 0.00967742, which seems sensible. 2021/12/16/06:00:48 PM Dropping snps with null values 2021/12/16/06:00:49 PM Metadata: 2021/12/16/06:00:49 PM Mean chi^2 = 2.838 2021/12/16/06:00:50 PM Lambda GC = 1.916 2021/12/16/06:00:50 PM Max chi^2 = 626.746 2021/12/16/06:00:50 PM 77637 Genome-wide significant SNPs (some may have been removed by filtering). 2021/12/16/06:00:50 PM Conversion finished at Thu Dec 16 18:00:50 2021 2021/12/16/06:00:50 PM Total time elapsed: 1.0m:29.08s 2021/12/16/06:01:01 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/06:01:01 PM Munging of Trait 2 complete. SNPs remaining: 7087911 2021/12/16/06:01:01 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><>

2021/12/16/06:01:13 PM Trait 2: Dropped 30 SNPs for duplicate values in the "snp_name" column 2021/12/16/06:01:18 PM Dropped 1039218 SNPs due to strand ambiguity, 5927870 SNPs remain in intersection after merging trait1 2021/12/16/06:01:39 PM Flipped the signs of of 3098396 SNPs to make them consistent with the effect allele orderings of the first trait. 2021/12/16/06:01:43 PM Dropped 0 SNPs due to strand ambiguity, 5670365 SNPs remain in intersection after merging trait2 2021/12/16/06:01:43 PM ... Merge of GWAS summary statistics complete. Number of SNPs: 5670365 2021/12/16/06:02:00 PM Using 5670365 SNPs to estimate Omega (0 SNPs excluded due to strand ambiguity) 2021/12/16/06:02:00 PM Estimating sigma.. 2021/12/16/06:03:02 PM Checking for positive definiteness .. 2021/12/16/06:03:02 PM Sigma hat: [[1.006 0.018] [0.018 1.078]] 2021/12/16/06:03:03 PM Beginning estimation of Omega ... 2021/12/16/06:03:03 PM Using GMM estimator of Omega .. 2021/12/16/06:03:03 PM Checking for positive definiteness .. 2021/12/16/06:03:03 PM Completed estimation of Omega ... 2021/12/16/06:03:03 PM Beginning MTAG calculations... 2021/12/16/06:03:10 PM ... Completed MTAG calculations. 2021/12/16/06:03:10 PM Writing Phenotype 1 to file ... 2021/12/16/06:04:06 PM Writing Phenotype 2 to file ... 2021/12/16/06:05:02 PM Summary of MTAG results:

Trait # SNPs used ... MTAG mean chi^2 GWAS equiv. (max) N 1 ...TRAIT_1 5670365 ... 1.180 233602 2 ...TRAIT_4 5670365 ... 2.679 747940

[2 rows x 7 columns]

Estimated Omega: [[7.761e-07 3.538e-07] [3.538e-07 2.452e-06]]

(Correlation): [[1. 0.256] [0.256 1. ]]

Estimated Sigma: [[1.006 0.018] [0.018 1.078]]

(Correlation): [[1. 0.017] [0.017 1. ]]

MTAG weight factors: (average across SNPs) [1.279 1.022]

2021/12/16/06:05:02 PM 2021/12/16/06:05:02 PM MTAG results saved to file. 2021/12/16/06:05:02 PM Beginning maxFDR calculations. Depending on the number of grid points specified, this might take some time... 2021/12/16/06:05:02 PM T=2 2021/12/16/06:05:02 PM Number of gridpoints to search: 200 2021/12/16/06:05:02 PM Performing grid search using 1 cores. 2021/12/16/06:05:02 PM Grid search: 10.0 percent finished for . Time: 0.001 min 2021/12/16/06:05:02 PM Grid search: 20.0 percent finished for . Time: 0.005 min 2021/12/16/06:05:02 PM Grid search: 30.0 percent finished for . Time: 0.007 min 2021/12/16/06:05:02 PM Grid search: 40.0 percent finished for . Time: 0.007 min 2021/12/16/06:05:02 PM Grid search: 50.0 percent finished for . Time: 0.008 min 2021/12/16/06:05:02 PM Grid search: 60.0 percent finished for . Time: 0.009 min 2021/12/16/06:05:03 PM Grid search: 70.0 percent finished for . Time: 0.010 min 2021/12/16/06:05:03 PM Grid search: 80.0 percent finished for . Time: 0.011 min 2021/12/16/06:05:03 PM Grid search: 90.0 percent finished for . Time: 0.012 min 2021/12/16/06:05:03 PM Grid search: 100.0 percent finished for . Time: 0.013 min 2021/12/16/06:05:03 PM Saved calculations of fdr over grid points in ./PAIRWISE_ANALYSIS_fdr_mat.txt 2021/12/16/06:05:03 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/06:05:03 PM grid point indices for max FDR for each trait: [99 39] 2021/12/16/06:05:03 PM Maximum FDR 2021/12/16/06:05:03 PM Max FDR of Trait 1: 0.0377427783459 at probs = [0.2 0.1 0.6 0.1] 2021/12/16/06:05:03 PM Max FDR of Trait 2: 6.90792700854e-06 at probs = [0. 0.5 0.3 0.2] 2021/12/16/06:05:03 PM <><><<>><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><><> 2021/12/16/06:05:03 PM Completed FDR calculations. 2021/12/16/06:05:03 PM MTAG complete. Time elapsed: 7.0m:57.301484108s

pjordab commented 2 years ago

Thank you very much!!!

paturley commented 2 years ago

It looks like several lines may be missing from your 2-trait log file. I don't see any lines about the maxFDR calculation.

On Fri, Dec 17, 2021 at 2:53 PM pjordab @.***> wrote:

Thank you very much!!!

— Reply to this email directly, view it on GitHub https://github.com/JonJala/mtag/issues/147#issuecomment-996994773, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFBUB5ISGBLVI2NC23GMNPDUROIKPANCNFSM5I27SLCQ . Triage notifications on the go with GitHub Mobile for iOS https://apps.apple.com/app/apple-store/id1477376905?ct=notification-email&mt=8&pt=524675 or Android https://play.google.com/store/apps/details?id=com.github.android&referrer=utm_campaign%3Dnotification-email%26utm_medium%3Demail%26utm_source%3Dgithub.

You are receiving this because you commented.Message ID: @.***>

pjordab commented 2 years ago

Oh sorry, I've updated it!

paturley commented 2 years ago

Alright. What I think is likely happening is that with more traits, there are more ways that the assumptions of MTAG can be violated. For example, with only 2 traits, each SNP can belong to one of 4 groups: associated with both, associated with trait 1 but not 2, associated with trait 2 but not 1, or not associated with either trait. MTAG's assmptions rely on every SNP belonging to either the first or last group. When you have 4 traits, there are 16 potential groups, but again, only two of them are consistent with MTAG's assumptions. Since this model is more flexible, it's possible to find more creative ways to inflate the FDR so the maxFDR will necessarily be larger when you have more traits.

pjordab commented 2 years ago

Thank you very much, this is really clarifying.

DonaldSandoz2000 commented 1 year ago

I I have the same problem, what should I do in this situation? Should I use only two traits for analysis, or what?

paturley commented 1 year ago

It depends on how comfortable you are with that high of an FDR. If you are OK with 40% of your hits potentially being false positives, then feel free to use them all. If you think a smaller FDR is important for your application, then I would use a smaller subset of traits such that the maxFDR is small enough.

On Wed, May 17, 2023 at 10:47 PM Donald Sandoz @.***> wrote:

I I have the same problem, what should I do in this situation? Should I use only two traits for analysis, or what?

— Reply to this email directly, view it on GitHub https://github.com/JonJala/mtag/issues/147#issuecomment-1552330952, or unsubscribe https://github.com/notifications/unsubscribe-auth/AFBUB5IQSY6G5VJ4GEOVOVLXGWEVFANCNFSM5I27SLCQ . You are receiving this because you commented.Message ID: @.***>