Atrial fibrillation - Heart Atrial Appendage
sheng Qian
2021-12-18
Last updated: 2022-02-08
Checks: 6 1
Knit directory: cTWAS_analysis/
This reproducible R Markdown analysis was created with workflowr (version 1.6.2). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
The command set.seed(20211220) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.
Using absolute paths to the files within your workflowr project makes it difficult for you and others to run your code on a different machine. Change the absolute path(s) below to the suggested relative path(s) to make your code more reproducible.
| absolute | relative |
|---|---|
| /project2/xinhe/shengqian/cTWAS/cTWAS_analysis/data/ | data |
| /project2/xinhe/shengqian/cTWAS/cTWAS_analysis/code/ctwas_config.R | code/ctwas_config.R |
Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.
The results in this page were generated with repository version ff0bb15. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.
Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .ipynb_checkpoints/
Ignored: analysis/figure/
Untracked files:
Untracked: analysis/BMI_Brain_Cortex.Rmd
Untracked: code/.ipynb_checkpoints/
Untracked: code/AF_Heart_Atrial_Appendage.err
Untracked: code/AF_Heart_Atrial_Appendage.out
Untracked: code/AF_Heart_Left_Ventricle.err
Untracked: code/AF_Heart_Left_Ventricle.out
Untracked: code/BMI_Brain_Cortex.err
Untracked: code/BMI_Brain_Cortex.out
Untracked: code/T2D_Adipose_Subcutaneous.err
Untracked: code/T2D_Adipose_Subcutaneous.out
Untracked: code/T2D_Adipose_Visceral_Omentum.err
Untracked: code/T2D_Adipose_Visceral_Omentum.out
Untracked: code/T2D_Liver.err
Untracked: code/T2D_Liver.out
Untracked: code/T2D_Pancreas.err
Untracked: code/T2D_Pancreas.out
Untracked: code/ctwas_config.R
Untracked: code/run_AF_analysis.sbatch
Untracked: code/run_AF_analysis.sh
Untracked: code/run_AF_ctwas_rss_LDR.R
Untracked: code/run_BMI_analysis.sbatch
Untracked: code/run_BMI_analysis.sh
Untracked: code/run_BMI_ctwas_rss_LDR.R
Untracked: code/run_T2D_analysis.sbatch
Untracked: code/run_T2D_analysis.sh
Untracked: code/run_T2D_ctwas_rss_LDR.R
Untracked: data/.ipynb_checkpoints/
Untracked: data/AF/
Untracked: data/BMI/
Untracked: data/T2D/
Untracked: data/UKB_460K.body_BMIz.sumstats
Untracked: data/Xue_et_al_T2D_META_Nat_Commun_2018
Untracked: data/gene_OMIM.txt
Untracked: data/gene_pip_0.8.txt
Untracked: data/mashr_Heart_Atrial_Appendage.db
Untracked: data/summary_known_genes_annotations.xlsx
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
These are the previous versions of the repository in which changes were made to the R Markdown (analysis/Atrial_Fibrillation_Heart_Atrial_Appendage.Rmd) and HTML (docs/Atrial_Fibrillation_Heart_Atrial_Appendage.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.
| File | Version | Author | Date | Message |
|---|---|---|---|---|
| Rmd | ff0bb15 | sq-96 | 2022-02-08 | update |
| html | 4106919 | sq-96 | 2021-12-24 | Build site. |
| html | f9aca3e | sq-96 | 2021-12-23 | Build site. |
| Rmd | d35d121 | sq-96 | 2021-12-23 | update analysis |
| html | 1f785cf | sq-96 | 2021-12-22 | Build site. |
| Rmd | 7784da6 | sq-96 | 2021-12-22 | Start my new project |
| html | 82c68fd | sq-96 | 2021-12-22 | Build site. |
| Rmd | 7e05bbc | sq-96 | 2021-12-22 | Start my new project |
| html | 9cf0e70 | sq-96 | 2021-12-22 | Build site. |
| Rmd | fbe71b7 | sq-96 | 2021-12-22 | Start my new project |
| html | 77de9fb | sq-96 | 2021-12-22 | Build site. |
| Rmd | c22d2e9 | sq-96 | 2021-12-22 | Start my new project |
| html | 6b46de7 | sq-96 | 2021-12-22 | Build site. |
| Rmd | 374327b | sq-96 | 2021-12-22 | Start my new project |
| html | e22d0c9 | sq-96 | 2021-12-21 | Build site. |
| html | a8dbae0 | sq-96 | 2021-12-21 | Build site. |
| Rmd | 37c4593 | sq-96 | 2021-12-21 | publish files |
| html | e080e7b | sq-96 | 2021-12-20 | Build site. |
| Rmd | 0b4a886 | sq-96 | 2021-12-20 | Start my new project |
| html | 9b1bb6e | sq-96 | 2021-12-20 | Build site. |
| html | 773dcb5 | sq-96 | 2021-12-20 | Build site. |
| html | 3052d49 | sq-96 | 2021-12-20 | Build site. |
| Rmd | 12069cf | sq-96 | 2021-12-20 | Publish the files |
Weight QC
qclist_all <- list()
qc_files <- paste0(results_dir, "/", list.files(results_dir, pattern="exprqc.Rd"))
for (i in 1:length(qc_files)){
load(qc_files[i])
chr <- unlist(strsplit(rev(unlist(strsplit(qc_files[i], "_")))[1], "[.]"))[1]
qclist_all[[chr]] <- cbind(do.call(rbind, lapply(qclist,unlist)), as.numeric(substring(chr,4)))
}
qclist_all <- data.frame(do.call(rbind, qclist_all))
colnames(qclist_all)[ncol(qclist_all)] <- "chr"
rm(qclist, wgtlist, z_gene_chr)
#number of imputed weights
nrow(qclist_all)[1] 10862#number of imputed weights by chromosome
table(qclist_all$chr)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
1103 776 625 427 501 624 518 377 399 423 648 611 203 368 356 495
17 18 19 20 21 22
676 177 840 318 126 271 #number of imputed weights without missing variants
sum(qclist_all$nmiss==0)[1] 7113#proportion of imputed weights without missing variants
mean(qclist_all$nmiss==0)[1] 0.6548518Load ctwas results
Check convergence of parameters
library(ggplot2)
library(cowplot)
********************************************************Note: As of version 1.0.0, cowplot does not change the default ggplot2 theme anymore. To recover the previous behavior, execute:
theme_set(theme_cowplot())********************************************************load(paste0(results_dir, "/", analysis_id, "_ctwas.s2.susieIrssres.Rd"))
df <- data.frame(niter = rep(1:ncol(group_prior_rec), 2),
value = c(group_prior_rec[1,], group_prior_rec[2,]),
group = rep(c("Gene", "SNP"), each = ncol(group_prior_rec)))
df$group <- as.factor(df$group)
df$value[df$group=="SNP"] <- df$value[df$group=="SNP"]*thin #adjust parameter to account for thin argument
p_pi <- ggplot(df, aes(x=niter, y=value, group=group)) +
geom_line(aes(color=group)) +
geom_point(aes(color=group)) +
xlab("Iteration") + ylab(bquote(pi)) +
ggtitle("Prior mean") +
theme_cowplot()
df <- data.frame(niter = rep(1:ncol(group_prior_var_rec), 2),
value = c(group_prior_var_rec[1,], group_prior_var_rec[2,]),
group = rep(c("Gene", "SNP"), each = ncol(group_prior_var_rec)))
df$group <- as.factor(df$group)
p_sigma2 <- ggplot(df, aes(x=niter, y=value, group=group)) +
geom_line(aes(color=group)) +
geom_point(aes(color=group)) +
xlab("Iteration") + ylab(bquote(sigma^2)) +
ggtitle("Prior variance") +
theme_cowplot()
plot_grid(p_pi, p_sigma2)
#estimated group prior
estimated_group_prior <- group_prior_rec[,ncol(group_prior_rec)]
names(estimated_group_prior) <- c("gene", "snp")
estimated_group_prior["snp"] <- estimated_group_prior["snp"]*thin #adjust parameter to account for thin argument
print(estimated_group_prior) gene snp
0.0208487453 0.0001927855 #estimated group prior variance
estimated_group_prior_var <- group_prior_var_rec[,ncol(group_prior_var_rec)]
names(estimated_group_prior_var) <- c("gene", "snp")
print(estimated_group_prior_var) gene snp
9.646456 8.679052 #report sample size
print(sample_size)[1] 1030836#report group size
group_size <- c(nrow(ctwas_gene_res), n_snps)
print(group_size)[1] 10862 6839050#estimated group PVE
estimated_group_pve <- estimated_group_prior_var*estimated_group_prior*group_size/sample_size #check PVE calculation
names(estimated_group_pve) <- c("gene", "snp")
print(estimated_group_pve) gene snp
0.00211918 0.01110076 #compare sum(PIP*mu2/sample_size) with above PVE calculation
c(sum(ctwas_gene_res$PVE),sum(ctwas_snp_res$PVE))[1] 0.01225055 0.12523867Genes with highest PIPs
#distribution of PIPs
hist(ctwas_gene_res$susie_pip, xlim=c(0,1), main="Distribution of Gene PIPs")
#genes with PIP>0.8 or 20 highest PIPs
head(ctwas_gene_res[order(-ctwas_gene_res$susie_pip),report_cols], max(sum(ctwas_gene_res$susie_pip>0.8), 20)) genename region_tag susie_pip mu2 PVE z
2293 CAV1 7_70 1.0000000 622.84815 6.042165e-04 15.567870
3275 PRRX1 1_84 0.9998808 119.34663 1.157627e-04 14.667578
4009 DEK 6_14 0.9925129 63.47784 6.111795e-05 -9.000000
3527 CCND2 12_4 0.9906174 27.01637 2.596231e-05 -5.283784
1310 PXN 12_75 0.9826498 29.47170 2.809405e-05 -5.328302
12836 RP11-325L7.2 5_82 0.9826020 1993.73789 1.900449e-03 12.356322
9257 AGAP5 10_49 0.9779202 48.91516 4.640420e-05 11.518590
3523 KLF12 13_36 0.9770829 26.00017 2.464439e-05 -5.072464
6621 AKAP6 14_8 0.9721549 76.81846 7.244551e-05 -9.197368
6914 JAM2 21_9 0.9640543 22.19279 2.075505e-05 4.563232
9639 DLEU1 13_21 0.9586468 23.59091 2.193885e-05 4.697095
11818 DPF3 14_34 0.9574312 33.35512 3.097993e-05 6.264960
10521 FAM43A 3_120 0.9569176 29.69898 2.756935e-05 -5.487179
2444 SEC23IP 10_74 0.9517809 22.38862 2.067163e-05 -4.565228
13075 LINC01629 14_36 0.9471613 31.88580 2.929757e-05 -5.695652
2138 AES 19_4 0.9403219 20.20438 1.843030e-05 4.182804
4658 POPDC3 6_70 0.9252511 25.19514 2.261449e-05 -4.758170
10290 NKX2-5 5_103 0.9228745 61.55975 5.511247e-05 -9.391892
7515 TNFSF13 17_7 0.9135188 34.26969 3.036953e-05 -5.883117
5185 GYPC 2_74 0.9068050 38.95864 3.427111e-05 -6.380531
8248 CMTM5 14_3 0.9067300 30.09803 2.647442e-05 -5.472727
10548 SCN10A 3_28 0.8867624 77.72022 6.685774e-05 -8.814286
13967 RP5-890E16.5 17_28 0.8660551 23.11200 1.941751e-05 -4.761194
712 SP100 2_135 0.8658723 18.73540 1.573719e-05 -3.671335
9012 MTSS1 8_82 0.8594108 20.87861 1.740655e-05 4.402634
8992 MURC 9_50 0.8478869 23.34736 1.920376e-05 4.911964
5223 PSMB7 9_64 0.8430652 25.30057 2.069197e-05 -4.820896
6114 STK11IP 2_130 0.8406381 18.96825 1.546845e-05 -3.868022
10416 PGP 16_2 0.8298586 28.51218 2.295329e-05 5.943820
9691 BOK 2_144 0.8255975 19.18259 1.536335e-05 3.910125
8420 MARS 12_36 0.8180336 17.81958 1.414097e-05 -3.366197
3088 GNB4 3_110 0.8097696 30.34628 2.383842e-05 -5.583333
num_eqtl
2293 3
3275 2
4009 1
3527 1
1310 1
12836 1
9257 2
3523 1
6621 1
6914 2
9639 1
11818 3
10521 1
2444 2
13075 1
2138 3
4658 1
10290 1
7515 1
5185 1
8248 1
10548 1
13967 1
712 2
9012 2
8992 2
5223 1
6114 2
10416 1
9691 3
8420 1
3088 1Genes with largest effect sizes
#plot PIP vs effect size
plot(ctwas_gene_res$susie_pip, ctwas_gene_res$mu2, xlab="PIP", ylab="mu^2", main="Gene PIPs vs Effect Size")
#genes with 20 largest effect sizes
head(ctwas_gene_res[order(-ctwas_gene_res$mu2),report_cols],20) genename region_tag susie_pip mu2 PVE z
3250 WIPF1 2_105 0.000000e+00 2852.2491 0.000000e+00 8.4415584
12836 RP11-325L7.2 5_82 9.826020e-01 1993.7379 1.900449e-03 12.3563218
1487 SIRT1 10_44 3.337242e-01 1299.2406 4.206179e-04 -5.0538462
7413 IL6R 1_75 5.862422e-12 1273.1807 7.240649e-15 -4.9784946
6444 HERC4 10_44 1.013394e-01 1220.9395 1.200281e-04 -5.3595598
11106 NACA 12_35 5.127731e-03 1179.0569 5.865032e-06 -6.2408092
960 BAZ2A 12_35 2.047295e-03 1162.5562 2.308898e-06 -5.9428571
7986 SLC35A1 6_59 0.000000e+00 804.8763 0.000000e+00 -5.0724638
2507 WNT3 17_27 3.707228e-05 797.4584 2.867925e-08 -4.3894024
2292 CAV2 7_70 0.000000e+00 689.1121 0.000000e+00 14.5349429
2293 CAV1 7_70 1.000000e+00 622.8482 6.042165e-04 15.5678701
4999 ORC3 6_59 0.000000e+00 597.6005 0.000000e+00 4.1666667
10616 ARL17A 17_27 7.713881e-06 488.6853 3.656896e-09 -0.5399498
3687 KDM3B 5_82 5.275983e-06 436.2575 2.232835e-09 -6.6287277
7728 GPR155 2_105 0.000000e+00 371.9331 0.000000e+00 -1.4546586
7276 ARHGAP27 17_27 5.875984e-06 189.4420 1.079860e-09 0.1576482
10739 MAPT 17_27 2.791224e-04 166.5127 4.508711e-08 3.8724621
7729 PMVK 1_76 9.693335e-05 161.1464 1.515319e-08 -12.1029412
7730 PBXIP1 1_76 9.335873e-05 156.1414 1.414111e-08 -11.8676471
341 FAM13B 5_82 1.826415e-01 153.9825 2.728233e-05 7.6009242
num_eqtl
3250 1
12836 1
1487 1
7413 1
6444 2
11106 2
960 1
7986 1
2507 2
2292 2
2293 3
4999 1
10616 2
3687 2
7728 2
7276 2
10739 2
7729 1
7730 1
341 2Genes with highest PVE
#genes with 20 highest pve
head(ctwas_gene_res[order(-ctwas_gene_res$PVE),report_cols],20) genename region_tag susie_pip mu2 PVE z
12836 RP11-325L7.2 5_82 0.9826020 1993.73789 1.900449e-03 12.356322
2293 CAV1 7_70 1.0000000 622.84815 6.042165e-04 15.567870
1487 SIRT1 10_44 0.3337242 1299.24064 4.206179e-04 -5.053846
6444 HERC4 10_44 0.1013394 1220.93950 1.200281e-04 -5.359560
3275 PRRX1 1_84 0.9998808 119.34663 1.157627e-04 14.667578
6621 AKAP6 14_8 0.9721549 76.81846 7.244551e-05 -9.197368
10548 SCN10A 3_28 0.8867624 77.72022 6.685774e-05 -8.814286
12013 ZSWIM8 10_49 0.7522384 87.25963 6.367652e-05 -11.216495
4009 DEK 6_14 0.9925129 63.47784 6.111795e-05 -9.000000
8298 SYNPO2L 10_49 0.6265676 99.72467 6.061512e-05 -11.945652
10290 NKX2-5 5_103 0.9228745 61.55975 5.511247e-05 -9.391892
9257 AGAP5 10_49 0.9779202 48.91516 4.640420e-05 11.518590
5185 GYPC 2_74 0.9068050 38.95864 3.427111e-05 -6.380531
3665 KCNJ5 11_80 0.5037532 67.73277 3.309993e-05 -8.748092
11818 DPF3 14_34 0.9574312 33.35512 3.097993e-05 6.264960
10515 C5orf47 5_104 0.7661297 41.61108 3.092586e-05 6.680556
7515 TNFSF13 17_7 0.9135188 34.26969 3.036953e-05 -5.883117
13075 LINC01629 14_36 0.9471613 31.88580 2.929757e-05 -5.695652
1310 PXN 12_75 0.9826498 29.47170 2.809405e-05 -5.328302
10521 FAM43A 3_120 0.9569176 29.69898 2.756935e-05 -5.487179
num_eqtl
12836 1
2293 3
1487 1
6444 2
3275 2
6621 1
10548 1
12013 1
4009 1
8298 1
10290 1
9257 2
5185 1
3665 1
11818 3
10515 1
7515 1
13075 1
1310 1
10521 1Genes with largest z scores
#genes with 20 largest z scores
head(ctwas_gene_res[order(-abs(ctwas_gene_res$z)),report_cols],20) genename region_tag susie_pip mu2 PVE z
2293 CAV1 7_70 1.000000e+00 622.84815 6.042165e-04 15.567870
3275 PRRX1 1_84 9.998808e-01 119.34663 1.157627e-04 14.667578
2292 CAV2 7_70 0.000000e+00 689.11215 0.000000e+00 14.534943
12836 RP11-325L7.2 5_82 9.826020e-01 1993.73789 1.900449e-03 12.356322
7729 PMVK 1_76 9.693335e-05 161.14638 1.515319e-08 -12.102941
8298 SYNPO2L 10_49 6.265676e-01 99.72467 6.061512e-05 -11.945652
7730 PBXIP1 1_76 9.335873e-05 156.14141 1.414111e-08 -11.867647
9805 MYOZ1 10_49 3.497772e-03 88.98387 3.019349e-07 11.759348
9257 AGAP5 10_49 9.779202e-01 48.91516 4.640420e-05 11.518590
12013 ZSWIM8 10_49 7.522384e-01 87.25963 6.367652e-05 -11.216495
6395 C9orf3 9_48 4.919856e-02 94.74228 4.521751e-06 10.671642
7407 ZBTB7B 1_76 9.636563e-05 122.09475 1.141378e-08 10.638889
13556 RP1-79C4.4 1_84 8.254209e-03 58.95140 4.720414e-07 9.586404
10290 NKX2-5 5_103 9.228745e-01 61.55975 5.511247e-05 -9.391892
9719 SEC24C 10_49 9.597200e-04 46.26092 4.306945e-08 9.271945
6621 AKAP6 14_8 9.721549e-01 76.81846 7.244551e-05 -9.197368
4009 DEK 6_14 9.925129e-01 63.47784 6.111795e-05 -9.000000
10548 SCN10A 3_28 8.867624e-01 77.72022 6.685774e-05 -8.814286
3665 KCNJ5 11_80 5.037532e-01 67.73277 3.309993e-05 -8.748092
7733 DCST2 1_76 2.613168e-04 89.76014 2.275418e-08 -8.718310
num_eqtl
2293 3
3275 2
2292 2
12836 1
7729 1
8298 1
7730 1
9805 2
9257 2
12013 1
6395 1
7407 1
13556 2
10290 1
9719 2
6621 1
4009 1
10548 1
3665 1
7733 1Comparing z scores and PIPs
#set nominal signifiance threshold for z scores
alpha <- 0.05
#bonferroni adjusted threshold for z scores
sig_thresh <- qnorm(1-(alpha/nrow(ctwas_gene_res)/2), lower=T)
#Q-Q plot for z scores
obs_z <- ctwas_gene_res$z[order(ctwas_gene_res$z)]
exp_z <- qnorm((1:nrow(ctwas_gene_res))/nrow(ctwas_gene_res))
plot(exp_z, obs_z, xlab="Expected z", ylab="Observed z", main="Gene z score Q-Q plot")
abline(a=0,b=1)
#plot z score vs PIP
plot(abs(ctwas_gene_res$z), ctwas_gene_res$susie_pip, xlab="abs(z)", ylab="PIP")
abline(v=sig_thresh, col="red", lty=2)
#proportion of significant z scores
mean(abs(ctwas_gene_res$z) > sig_thresh)[1] 0.008654023#genes with most significant z scores
head(ctwas_gene_res[order(-abs(ctwas_gene_res$z)),report_cols],20) genename region_tag susie_pip mu2 PVE z
2293 CAV1 7_70 1.000000e+00 622.84815 6.042165e-04 15.567870
3275 PRRX1 1_84 9.998808e-01 119.34663 1.157627e-04 14.667578
2292 CAV2 7_70 0.000000e+00 689.11215 0.000000e+00 14.534943
12836 RP11-325L7.2 5_82 9.826020e-01 1993.73789 1.900449e-03 12.356322
7729 PMVK 1_76 9.693335e-05 161.14638 1.515319e-08 -12.102941
8298 SYNPO2L 10_49 6.265676e-01 99.72467 6.061512e-05 -11.945652
7730 PBXIP1 1_76 9.335873e-05 156.14141 1.414111e-08 -11.867647
9805 MYOZ1 10_49 3.497772e-03 88.98387 3.019349e-07 11.759348
9257 AGAP5 10_49 9.779202e-01 48.91516 4.640420e-05 11.518590
12013 ZSWIM8 10_49 7.522384e-01 87.25963 6.367652e-05 -11.216495
6395 C9orf3 9_48 4.919856e-02 94.74228 4.521751e-06 10.671642
7407 ZBTB7B 1_76 9.636563e-05 122.09475 1.141378e-08 10.638889
13556 RP1-79C4.4 1_84 8.254209e-03 58.95140 4.720414e-07 9.586404
10290 NKX2-5 5_103 9.228745e-01 61.55975 5.511247e-05 -9.391892
9719 SEC24C 10_49 9.597200e-04 46.26092 4.306945e-08 9.271945
6621 AKAP6 14_8 9.721549e-01 76.81846 7.244551e-05 -9.197368
4009 DEK 6_14 9.925129e-01 63.47784 6.111795e-05 -9.000000
10548 SCN10A 3_28 8.867624e-01 77.72022 6.685774e-05 -8.814286
3665 KCNJ5 11_80 5.037532e-01 67.73277 3.309993e-05 -8.748092
7733 DCST2 1_76 2.613168e-04 89.76014 2.275418e-08 -8.718310
num_eqtl
2293 3
3275 2
2292 2
12836 1
7729 1
8298 1
7730 1
9805 2
9257 2
12013 1
6395 1
7407 1
13556 2
10290 1
9719 2
6621 1
4009 1
10548 1
3665 1
7733 1SNPs with highest PIPs
#snps with PIP>0.8 or 20 highest PIPs
head(ctwas_snp_res[order(-ctwas_snp_res$susie_pip),report_cols_snps],
max(sum(ctwas_snp_res$susie_pip>0.8), 20)) id region_tag susie_pip mu2 PVE z
28439 rs12404927 1_75 1.0000000 1437.87698 1.394865e-03 -1.5517241
28440 rs7536152 1_75 1.0000000 1473.01903 1.428956e-03 -6.4626866
28533 rs34515871 1_76 1.0000000 299.28422 2.903316e-04 17.0000000
32410 rs12142529 1_83 1.0000000 145.61337 1.412576e-04 12.0396040
32435 rs112797273 1_83 1.0000000 101.67077 9.862944e-05 9.4903846
183533 rs1906615 4_72 1.0000000 1756.73093 1.704181e-03 45.1604938
183538 rs7440714 4_72 1.0000000 992.10110 9.624238e-04 -9.2183406
272612 rs111990232 6_59 1.0000000 1710.90898 1.659730e-03 0.5392670
414968 rs76106073 10_44 1.0000000 1273.72173 1.235620e-03 0.6572238
550561 rs140798119 15_35 1.0000000 285.26104 2.767279e-04 3.9950739
550564 rs74022964 15_35 1.0000000 342.92230 3.326643e-04 12.5777778
681435 rs4972703 2_105 1.0000000 4299.78452 4.171163e-03 -0.7842670
739005 rs9770220 7_70 1.0000000 835.76918 8.107683e-04 -5.2814815
788780 rs61937778 12_35 1.0000000 1668.54978 1.618637e-03 -1.3056995
28537 rs2878412 1_76 1.0000000 105.57502 1.024169e-04 5.4250000
586354 rs62056842 17_27 1.0000000 1130.21900 1.096410e-03 1.1269036
28535 rs11576820 1_76 1.0000000 133.71597 1.297160e-04 8.0243902
422177 rs60469668 10_66 1.0000000 172.10824 1.669599e-04 16.9647059
110779 rs9852222 3_9 1.0000000 83.84308 8.133503e-05 8.3823529
339255 rs373983 8_21 0.9999998 67.64329 6.561982e-05 8.5588235
298378 rs12112152 7_15 0.9999995 45.13406 4.378392e-05 6.9358974
422175 rs7094488 10_66 0.9999992 110.01432 1.067233e-04 14.1060606
28540 rs10908445 1_76 0.9999991 193.56990 1.877794e-04 -13.1470588
570907 rs6499606 16_39 0.9999990 103.94144 1.008321e-04 12.6376812
826958 rs8005417 14_9 0.9999978 427.56249 4.147716e-04 0.9863946
852634 rs140185678 16_2 0.9999910 45.86064 4.448839e-05 7.6100917
97513 rs1975584 2_118 0.9999847 71.76048 6.961280e-05 1.2311828
580602 rs72811292 17_11 0.9999789 40.04445 3.884576e-05 -6.5233645
279525 rs11756438 6_79 0.9999723 47.69440 4.626641e-05 -8.2089552
117531 rs116202356 3_27 0.9999305 36.10604 3.502355e-05 6.0813559
176218 rs1458038 4_54 0.9998791 34.99562 3.394467e-05 6.0277778
121102 rs12330500 3_40 0.9997290 519.11865 5.034535e-04 0.5641026
396844 rs1886296 9_73 0.9995279 34.54923 3.349992e-05 5.9863014
39581 rs6427989 1_102 0.9993189 36.39724 3.528442e-05 -4.4142857
653873 rs464901 22_4 0.9992684 43.84562 4.250292e-05 -7.0555556
708945 rs574293775 5_82 0.9987918 1853.69252 1.796069e-03 6.2411765
39582 rs12353975 1_102 0.9982647 35.87479 3.474126e-05 -4.3000000
184120 rs7700110 4_73 0.9979235 45.19125 4.374839e-05 7.0666667
232199 rs199992924 5_68 0.9969497 1340.40151 1.296339e-03 -0.6309524
655632 rs133902 22_7 0.9963501 31.05883 3.001978e-05 6.1617647
232019 rs338623 5_68 0.9953592 67.12738 6.481716e-05 -7.9855072
414975 rs12360521 10_44 0.9953399 1283.59990 1.239400e-03 5.3740458
696421 rs7374540 3_28 0.9939021 45.86079 4.421764e-05 4.7794118
570915 rs876727 16_39 0.9926586 48.82838 4.702000e-05 -10.0000000
275680 rs9496567 6_67 0.9918710 26.51232 2.551017e-05 5.1125000
586570 rs75230966 17_27 0.9917637 49.43748 4.756362e-05 6.0582524
183553 rs4631108 4_72 0.9915393 306.03072 2.943645e-04 -17.6376812
40258 rs4951023 1_104 0.9906227 31.13453 2.991996e-05 5.5970149
238444 rs6894302 5_84 0.9877503 34.35180 3.291601e-05 -7.2933333
666039 rs2885697 1_25 0.9874950 44.28194 4.242013e-05 -6.2714286
570920 rs60602157 16_39 0.9854801 38.72262 3.701885e-05 10.4415584
339226 rs403894 8_19 0.9793925 56.27052 5.346236e-05 6.8529412
650809 rs7282237 21_16 0.9782270 25.00237 2.372636e-05 -4.7605634
529120 rs2738413 14_29 0.9782089 112.52732 1.067825e-04 -11.6119403
570912 rs4788691 16_39 0.9747904 40.13308 3.795108e-05 -1.6900000
329588 rs35760656 7_94 0.9746629 35.32804 3.340291e-05 -6.3536585
470281 rs17380837 12_18 0.9737999 44.42786 4.196967e-05 -6.9583333
366399 rs7460121 8_88 0.9725824 24.19453 2.282727e-05 4.7500000
2727 rs284278 1_7 0.9713320 35.00958 3.298869e-05 -6.0857143
826964 rs8011559 14_9 0.9693907 435.54854 4.095867e-04 5.6913580
283467 rs958747 6_89 0.9688123 24.37662 2.290992e-05 -4.9710145
556130 rs12898337 15_48 0.9643924 35.44339 3.315885e-05 -6.0579710
279482 rs77435894 6_78 0.9642414 28.80796 2.694690e-05 6.1111111
549931 rs745636 15_33 0.9571831 25.61925 2.378876e-05 -5.0740741
25783 rs4839174 1_69 0.9569309 33.97079 3.153528e-05 5.9459459
183505 rs1823291 4_72 0.9568146 286.06204 2.655207e-04 -19.4520548
604563 rs17794590 18_24 0.9550117 24.59348 2.278448e-05 -4.8242424
100283 rs35880620 2_125 0.9544571 39.77143 3.682460e-05 -6.2168675
58743 rs7578482 2_16 0.9542037 32.14647 2.975670e-05 5.7435897
738958 rs1633714 7_70 0.9499686 778.22440 7.171740e-04 11.7536232
97488 rs11889306 2_118 0.9487414 45.59873 4.196730e-05 6.7391304
826887 rs73241997 14_9 0.9472521 276.10117 2.537139e-04 7.8817204
681445 rs1367220 2_105 0.9467808 4303.78606 3.952852e-03 7.2207792
784225 rs2291437 12_17 0.9462388 66.11326 6.068757e-05 9.1826923
32457 rs7522387 1_83 0.9460866 59.86066 5.493926e-05 6.6375000
696586 rs7373065 3_28 0.9355464 65.79596 5.971384e-05 -8.0637450
28495 rs906280 1_75 0.9353945 59.57170 5.405616e-05 7.9111111
490521 rs12425471 12_69 0.9303264 41.49244 3.744680e-05 -5.2500000
614935 rs7256735 19_3 0.9241700 23.45925 2.103180e-05 4.4952381
279432 rs89107 6_78 0.9189049 47.36791 4.222457e-05 -9.5606061
372071 rs1594768 9_9 0.9139999 24.21822 2.147330e-05 -4.6911765
5627 rs10917072 1_15 0.9043759 33.53922 2.942472e-05 5.8235294
217602 rs114414434 5_31 0.8994692 26.16276 2.282866e-05 -4.9718310
232192 rs4073838 5_68 0.8923424 1368.92743 1.185011e-03 -5.9545455
272609 rs371814 6_59 0.8922509 1722.76723 1.491159e-03 6.1911765
252160 rs73724866 6_13 0.8855608 94.00847 8.075990e-05 -10.4141414
479761 rs776211 12_43 0.8831880 28.12821 2.409937e-05 -6.3424658
550553 rs519946 15_34 0.8819650 34.47055 2.949239e-05 5.9242424
494561 rs6560886 12_82 0.8817026 28.94415 2.475673e-05 5.6666667
28531 rs12128882 1_76 0.8779576 32.56646 2.773668e-05 1.0731707
68039 rs243080 2_40 0.8749115 24.16683 2.051135e-05 -4.5522388
467254 rs12821447 12_12 0.8716947 24.08112 2.036346e-05 4.5128205
615414 rs73919353 19_5 0.8614285 28.38347 2.371894e-05 5.4076433
453904 rs565449 11_55 0.8575434 28.38859 2.361622e-05 -5.2535211
123427 rs1091584 3_45 0.8532454 30.27502 2.505929e-05 5.5882353
244628 rs62377226 5_100 0.8422335 33.87890 2.768039e-05 5.8493151
113948 rs73032373 3_18 0.8375231 33.29814 2.705373e-05 -5.9178082
369034 rs72693377 8_94 0.8312985 24.72583 1.993968e-05 4.2762431
659021 rs11705586 22_15 0.8290134 26.32683 2.117242e-05 4.8988764
121105 rs6924 3_40 0.8277296 522.73041 4.197364e-04 4.2133333
304278 rs145593380 7_26 0.8170209 29.74581 2.357596e-05 5.3333333
495691 rs7321083 13_3 0.8149736 24.65878 1.949510e-05 4.8292683SNPs with largest effect sizes
#plot PIP vs effect size
#plot(ctwas_snp_res$susie_pip, ctwas_snp_res$mu2, xlab="PIP", ylab="mu^2", main="SNP PIPs vs Effect Size")
#SNPs with 50 largest effect sizes
head(ctwas_snp_res[order(-ctwas_snp_res$mu2),report_cols_snps],50) id region_tag susie_pip mu2 PVE z
681445 rs1367220 2_105 9.467808e-01 4303.786 3.952852e-03 7.220779
681435 rs4972703 2_105 1.000000e+00 4299.785 4.171163e-03 -0.784267
681446 rs1367219 2_105 2.069799e-01 4294.507 8.622874e-04 7.103896
681441 rs10168156 2_105 3.312313e-01 4292.416 1.379252e-03 7.233766
681439 rs6713018 2_105 9.687858e-02 4291.538 4.033213e-04 7.168831
681430 rs1864453 2_105 1.876497e-01 4289.961 7.809293e-04 7.207792
681422 rs2033315 2_105 1.145049e-01 4287.502 4.762543e-04 7.064935
681433 rs6707162 2_105 3.496928e-02 4279.920 1.451887e-04 7.090909
681434 rs6735680 2_105 2.455791e-02 4279.461 1.019509e-04 7.077922
681423 rs2033314 2_105 2.432678e-02 4278.387 1.009660e-04 7.155844
681425 rs28485554 2_105 8.747105e-03 4277.770 3.629879e-05 7.000000
681456 rs10803884 2_105 7.083440e-03 4217.525 2.898093e-05 7.285714
681464 rs6738901 2_105 7.743316e-03 4201.159 3.155778e-05 -7.272727
681449 rs12466643 2_105 3.135672e-04 4196.140 1.276412e-06 7.272727
681480 rs35368253 2_105 1.998443e-04 4178.444 8.100592e-07 -7.256410
681466 rs10197521 2_105 1.763290e-11 4171.648 7.135785e-14 -6.547619
681468 rs34661753 2_105 3.194389e-08 4119.395 1.276532e-10 -7.076923
681419 rs7590328 2_105 2.272449e-11 4103.428 9.045890e-14 6.784810
681417 rs6433497 2_105 0.000000e+00 3465.041 0.000000e+00 6.345679
681418 rs2115874 2_105 0.000000e+00 3462.201 0.000000e+00 6.419753
681410 rs1430185 2_105 0.000000e+00 3387.075 0.000000e+00 6.185185
681397 rs1991601 2_105 0.000000e+00 3363.058 0.000000e+00 5.950617
681395 rs6759870 2_105 0.000000e+00 2978.209 0.000000e+00 5.987342
681486 rs35215597 2_105 0.000000e+00 2793.497 0.000000e+00 -8.441558
681506 rs56181519 2_105 0.000000e+00 2703.083 0.000000e+00 -8.597403
681308 rs13032076 2_105 0.000000e+00 2402.122 0.000000e+00 -5.250000
681408 rs13024657 2_105 0.000000e+00 2366.705 0.000000e+00 -5.258065
681509 rs2358891 2_105 0.000000e+00 2265.703 0.000000e+00 -7.500000
681526 rs35808589 2_105 0.000000e+00 2030.164 0.000000e+00 -5.229885
681558 rs1376875 2_105 0.000000e+00 2019.770 0.000000e+00 -5.000000
681541 rs35444726 2_105 0.000000e+00 2016.551 0.000000e+00 -5.172414
681550 rs13420492 2_105 0.000000e+00 1989.902 0.000000e+00 -4.908046
681552 rs71417497 2_105 0.000000e+00 1980.143 0.000000e+00 -4.772727
708976 rs2040862 5_82 1.932759e-01 1961.536 3.677770e-04 12.459770
708860 rs17171711 5_82 5.061870e-02 1958.969 9.619419e-05 12.482759
708884 rs13355516 5_82 4.398123e-03 1956.728 8.348495e-06 12.379310
708835 rs9327807 5_82 2.294241e-03 1953.970 4.348779e-06 12.367816
708869 rs77915370 5_82 3.425950e-04 1952.684 6.489682e-07 12.356322
708740 rs78081438 5_82 6.691465e-04 1951.053 1.266487e-06 12.390805
708691 rs73300168 5_82 3.565820e-04 1950.786 6.748068e-07 12.298851
708817 rs11959181 5_82 9.541911e-05 1949.634 1.804674e-07 12.241379
708671 rs10076361 5_82 2.676249e-04 1949.238 5.060598e-07 12.298851
708836 rs113331339 5_82 8.891833e-05 1949.202 1.681352e-07 12.241379
708789 rs73299268 5_82 8.169862e-05 1949.085 1.544742e-07 12.241379
708892 rs148378888 5_82 1.992837e-04 1948.655 3.767186e-07 12.344828
708714 rs73299210 5_82 5.888915e-05 1946.246 1.111843e-07 12.287356
708717 rs73299219 5_82 6.204555e-05 1945.825 1.171183e-07 12.287356
708737 rs13362264 5_82 2.545517e-05 1945.354 4.803802e-08 12.264368
681335 rs2303891 2_105 0.000000e+00 1921.675 0.000000e+00 -4.185185
708899 rs10479176 5_82 6.731608e-01 1866.560 1.218908e-03 6.931818SNPs with highest PVE
#SNPs with 50 highest pve
head(ctwas_snp_res[order(-ctwas_snp_res$PVE),report_cols_snps],50) id region_tag susie_pip mu2 PVE z
681435 rs4972703 2_105 1.00000000 4299.7845 0.0041711625 -0.7842670
681445 rs1367220 2_105 0.94678083 4303.7861 0.0039528520 7.2207792
708945 rs574293775 5_82 0.99879179 1853.6925 0.0017960693 6.2411765
183533 rs1906615 4_72 1.00000000 1756.7309 0.0017041808 45.1604938
272612 rs111990232 6_59 1.00000000 1710.9090 0.0016597296 0.5392670
788780 rs61937778 12_35 1.00000000 1668.5498 0.0016186375 -1.3056995
272609 rs371814 6_59 0.89225093 1722.7672 0.0014911593 6.1911765
28440 rs7536152 1_75 1.00000000 1473.0190 0.0014289557 -6.4626866
28439 rs12404927 1_75 1.00000000 1437.8770 0.0013948649 -1.5517241
681441 rs10168156 2_105 0.33123129 4292.4155 0.0013792517 7.2337662
232199 rs199992924 5_68 0.99694972 1340.4015 0.0012963390 -0.6309524
414975 rs12360521 10_44 0.99533986 1283.5999 0.0012394000 5.3740458
414968 rs76106073 10_44 1.00000000 1273.7217 0.0012356201 0.6572238
708899 rs10479176 5_82 0.67316075 1866.5597 0.0012189085 6.9318182
232192 rs4073838 5_68 0.89234244 1368.9274 0.0011850110 -5.9545455
586354 rs62056842 17_27 1.00000000 1130.2190 0.0010964101 1.1269036
788720 rs2860482 12_35 0.64493512 1697.2111 0.0010618479 -7.1052632
788775 rs7313074 12_35 0.58220672 1705.1065 0.0009630285 -6.9868421
183538 rs7440714 4_72 1.00000000 992.1011 0.0009624238 -9.2183406
681446 rs1367219 2_105 0.20697994 4294.5073 0.0008622874 7.1038961
739005 rs9770220 7_70 1.00000000 835.7692 0.0008107683 -5.2814815
681430 rs1864453 2_105 0.18764974 4289.9609 0.0007809293 7.2077922
738958 rs1633714 7_70 0.94996864 778.2244 0.0007171740 11.7536232
788778 rs4759256 12_35 0.41979796 1704.6254 0.0006941921 -6.9736842
272613 rs384318 6_59 0.37976257 1721.2083 0.0006340975 6.1029412
738995 rs10255816 7_70 0.78914390 709.0714 0.0005428210 -9.7123288
121102 rs12330500 3_40 0.99972902 519.1187 0.0005034535 0.5641026
681422 rs2033315 2_105 0.11450492 4287.5023 0.0004762543 7.0649351
272614 rs1145714 6_59 0.26478699 1717.6032 0.0004411943 6.1029412
121105 rs6924 3_40 0.82772961 522.7304 0.0004197364 4.2133333
826958 rs8005417 14_9 0.99999777 427.5625 0.0004147716 0.9863946
826964 rs8011559 14_9 0.96939072 435.5485 0.0004095867 5.6913580
681439 rs6713018 2_105 0.09687858 4291.5379 0.0004033213 7.1688312
232191 rs4235764 5_68 0.29136909 1367.2078 0.0003864457 -5.8787879
708976 rs2040862 5_82 0.19327595 1961.5362 0.0003677770 12.4597701
788715 rs7978685 12_35 0.21731446 1694.8625 0.0003573004 -7.0657895
550564 rs74022964 15_35 1.00000000 342.9223 0.0003326643 12.5777778
232194 rs4235768 5_68 0.24121476 1367.3620 0.0003199616 -5.8030303
272571 rs9444476 6_59 0.19008883 1698.6999 0.0003132447 -6.3088235
183553 rs4631108 4_72 0.99153932 306.0307 0.0002943645 -17.6376812
28533 rs34515871 1_76 1.00000000 299.2842 0.0002903316 17.0000000
272596 rs9444488 6_59 0.16974071 1703.6881 0.0002805347 -6.2647059
550561 rs140798119 15_35 1.00000000 285.2610 0.0002767279 3.9950739
183505 rs1823291 4_72 0.95681461 286.0620 0.0002655207 -19.4520548
826887 rs73241997 14_9 0.94725212 276.1012 0.0002537139 7.8817204
708944 rs191830974 5_82 0.12101272 1865.7774 0.0002190288 6.9206349
121104 rs1554125 3_40 0.38957721 520.5125 0.0001967139 4.1333333
28540 rs10908445 1_76 0.99999905 193.5699 0.0001877794 -13.1470588
272583 rs9444481 6_59 0.10849613 1701.3375 0.0001790668 -6.2500000
422177 rs60469668 10_66 1.00000000 172.1082 0.0001669599 16.9647059SNPs with largest z scores
#histogram of (abs) SNP z scores
hist(abs(ctwas_snp_res$z))
#SNPs with 50 largest z scores
head(ctwas_snp_res[order(-abs(ctwas_snp_res$z)),report_cols_snps],50) id region_tag susie_pip mu2 PVE z
183533 rs1906615 4_72 1.000000e+00 1756.7309 1.704181e-03 45.16049
183534 rs75725917 4_72 0.000000e+00 588.7458 0.000000e+00 42.33663
183532 rs1906611 4_72 0.000000e+00 526.9126 0.000000e+00 40.97143
183531 rs28521134 4_72 0.000000e+00 512.6248 0.000000e+00 40.22772
183527 rs10019689 4_72 0.000000e+00 508.5039 0.000000e+00 40.11000
183528 rs76013973 4_72 0.000000e+00 502.3963 0.000000e+00 40.04902
183530 rs12639820 4_72 0.000000e+00 499.7947 0.000000e+00 40.03922
183526 rs12647393 4_72 0.000000e+00 500.8094 0.000000e+00 39.98020
183529 rs74496596 4_72 0.000000e+00 484.7331 0.000000e+00 39.75728
183525 rs12647316 4_72 0.000000e+00 476.7021 0.000000e+00 39.55882
183524 rs4529121 4_72 0.000000e+00 480.5177 0.000000e+00 39.55340
183518 rs12650829 4_72 0.000000e+00 450.0533 0.000000e+00 31.54444
183541 rs3866831 4_72 0.000000e+00 919.4634 0.000000e+00 -30.10000
183540 rs6533530 4_72 0.000000e+00 904.8879 0.000000e+00 -29.91429
183523 rs12644107 4_72 0.000000e+00 153.0712 0.000000e+00 25.98851
183521 rs2723318 4_72 0.000000e+00 365.6985 0.000000e+00 -25.45833
183516 rs2218698 4_72 0.000000e+00 337.5178 0.000000e+00 -25.09589
183519 rs2197814 4_72 0.000000e+00 339.8970 0.000000e+00 24.59722
183517 rs1448799 4_72 0.000000e+00 339.2379 0.000000e+00 24.56944
183515 rs112927894 4_72 0.000000e+00 320.7713 0.000000e+00 24.50685
183505 rs1823291 4_72 9.568146e-01 286.0620 2.655207e-04 -19.45205
183509 rs2723296 4_72 4.318269e-02 278.0630 1.164832e-05 -19.02740
183514 rs11724067 4_72 0.000000e+00 321.1887 0.000000e+00 -18.94059
183508 rs2044674 4_72 2.641186e-06 253.4848 6.494734e-10 18.58108
183553 rs4631108 4_72 9.915393e-01 306.0307 2.943645e-04 -17.63768
183552 rs1906613 4_72 8.460684e-03 295.1509 2.422479e-06 17.34783
28533 rs34515871 1_76 1.000000e+00 299.2842 2.903316e-04 17.00000
422177 rs60469668 10_66 1.000000e+00 172.1082 1.669599e-04 16.96471
183513 rs13111704 4_72 0.000000e+00 393.8288 0.000000e+00 -16.91743
28536 rs12058931 1_76 1.643627e-04 221.4339 3.530676e-08 16.54054
183537 rs4124159 4_72 0.000000e+00 1159.1777 0.000000e+00 -16.14414
183535 rs1906606 4_72 0.000000e+00 1154.4323 0.000000e+00 -16.02703
183542 rs4032974 4_72 0.000000e+00 1146.0129 0.000000e+00 -15.86726
739346 rs11773845 7_70 0.000000e+00 555.7324 0.000000e+00 15.73134
183562 rs17513625 4_72 0.000000e+00 176.0291 0.000000e+00 15.71130
183539 rs10006881 4_72 0.000000e+00 1136.3170 0.000000e+00 -15.68750
739366 rs1997571 7_70 0.000000e+00 553.9040 0.000000e+00 15.68657
739329 rs3807989 7_70 0.000000e+00 551.5677 0.000000e+00 15.62687
422173 rs12572965 10_66 2.273800e-04 132.7631 2.928466e-08 15.27059
422172 rs56965730 10_66 2.300323e-04 131.7918 2.940950e-08 15.21176
739367 rs1997572 7_70 0.000000e+00 525.1335 0.000000e+00 14.50685
739223 rs2270188 7_70 0.000000e+00 711.6787 0.000000e+00 14.43939
739224 rs2270189 7_70 0.000000e+00 711.6066 0.000000e+00 14.43939
739272 rs2109514 7_70 0.000000e+00 687.7806 0.000000e+00 14.42424
739276 rs55883210 7_70 0.000000e+00 687.5295 0.000000e+00 14.42424
739239 rs10271007 7_70 0.000000e+00 707.4203 0.000000e+00 14.39394
739256 rs6466579 7_70 0.000000e+00 703.1231 0.000000e+00 14.33333
739267 rs7795510 7_70 0.000000e+00 696.7170 0.000000e+00 14.22727
183507 rs12642151 4_72 2.282285e-12 237.6299 5.261159e-16 14.14286
422175 rs7094488 10_66 9.999992e-01 110.0143 1.067233e-04 14.10606Gene set enrichment for genes with PIP>0.8
#GO enrichment analysis
library(enrichR)Welcome to enrichR
Checking connection ... Enrichr ... Connection is Live!
FlyEnrichr ... Connection is available!
WormEnrichr ... Connection is available!
YeastEnrichr ... Connection is available!
FishEnrichr ... Connection is available!dbs <- c("GO_Biological_Process_2021", "GO_Cellular_Component_2021", "GO_Molecular_Function_2021")
genes <- ctwas_gene_res$genename[ctwas_gene_res$susie_pip>0.5]
#number of genes for gene set enrichment
length(genes)[1] 68if (length(genes)>0){
GO_enrichment <- enrichr(genes, dbs)
for (db in dbs){
print(db)
df <- GO_enrichment[[db]]
print(plotEnrich(GO_enrichment[[db]]))
df <- df[df$Adjusted.P.value<0.05,c("Term", "Overlap", "Adjusted.P.value", "Genes")]
print(df)
}
#DisGeNET enrichment
# devtools::install_bitbucket("ibi_group/disgenet2r")
library(disgenet2r)
disgenet_api_key <- get_disgenet_api_key(
email = "wesleycrouse@gmail.com",
password = "uchicago1" )
Sys.setenv(DISGENET_API_KEY= disgenet_api_key)
res_enrich <-disease_enrichment(entities=genes, vocabulary = "HGNC",
database = "CURATED" )
df <- res_enrich@qresult[1:10, c("Description", "FDR", "Ratio", "BgRatio")]
print(df)
#WebGestalt enrichment
library(WebGestaltR)
background <- ctwas_gene_res$genename
#listGeneSet()
databases <- c("pathway_KEGG", "disease_GLAD4U", "disease_OMIM")
enrichResult <- WebGestaltR(enrichMethod="ORA", organism="hsapiens",
interestGene=genes, referenceGene=background,
enrichDatabase=databases, interestGeneType="genesymbol",
referenceGeneType="genesymbol", isOutput=F)
print(enrichResult[,c("description", "size", "overlap", "FDR", "database", "userId")])
}Uploading data to Enrichr... Done.
Querying GO_Biological_Process_2021... Done.
Querying GO_Cellular_Component_2021... Done.
Querying GO_Molecular_Function_2021... Done.
Parsing results... Done.
[1] "GO_Biological_Process_2021"
| Version | Author | Date |
|---|---|---|
| 77de9fb | sq-96 | 2021-12-22 |
[1] Term Overlap Adjusted.P.value Genes
<0 rows> (or 0-length row.names)
[1] "GO_Cellular_Component_2021"
Term Overlap Adjusted.P.value
1 plasma membrane raft (GO:0044853) 4/82 0.01680691
Genes
1 SCARB1;STIM1;CAV1;AKAP6
[1] "GO_Molecular_Function_2021"[1] Term Overlap Adjusted.P.value Genes
<0 rows> (or 0-length row.names)AHSA2 gene(s) from the input list not found in DisGeNET CURATEDZBTB38 gene(s) from the input list not found in DisGeNET CURATEDSEC23IP gene(s) from the input list not found in DisGeNET CURATEDEXOC8 gene(s) from the input list not found in DisGeNET CURATEDCMTM5 gene(s) from the input list not found in DisGeNET CURATEDPWP1 gene(s) from the input list not found in DisGeNET CURATEDRAB1A gene(s) from the input list not found in DisGeNET CURATEDNFATC2IP gene(s) from the input list not found in DisGeNET CURATEDZSWIM8 gene(s) from the input list not found in DisGeNET CURATEDBOK gene(s) from the input list not found in DisGeNET CURATEDMURC gene(s) from the input list not found in DisGeNET CURATEDECE2 gene(s) from the input list not found in DisGeNET CURATEDDLEU1 gene(s) from the input list not found in DisGeNET CURATEDGLTSCR1 gene(s) from the input list not found in DisGeNET CURATEDPAFAH1B2 gene(s) from the input list not found in DisGeNET CURATEDJAM2 gene(s) from the input list not found in DisGeNET CURATEDGMCL1 gene(s) from the input list not found in DisGeNET CURATEDLINC01629 gene(s) from the input list not found in DisGeNET CURATEDGIT2 gene(s) from the input list not found in DisGeNET CURATEDAGAP5 gene(s) from the input list not found in DisGeNET CURATEDMARS gene(s) from the input list not found in DisGeNET CURATEDC5orf47 gene(s) from the input list not found in DisGeNET CURATEDSP100 gene(s) from the input list not found in DisGeNET CURATEDAES gene(s) from the input list not found in DisGeNET CURATEDC9orf43 gene(s) from the input list not found in DisGeNET CURATEDRP11-325L7.2 gene(s) from the input list not found in DisGeNET CURATEDNCLN gene(s) from the input list not found in DisGeNET CURATEDFAM57A gene(s) from the input list not found in DisGeNET CURATEDPOPDC3 gene(s) from the input list not found in DisGeNET CURATEDSTK11IP gene(s) from the input list not found in DisGeNET CURATEDNBL1 gene(s) from the input list not found in DisGeNET CURATEDLIN54 gene(s) from the input list not found in DisGeNET CURATEDRP5-890E16.5 gene(s) from the input list not found in DisGeNET CURATED Description FDR Ratio BgRatio
135 Paroxysmal atrial fibrillation 3.461631e-13 13/35 156/9703
257 Persistent atrial fibrillation 3.461631e-13 13/35 156/9703
279 familial atrial fibrillation 3.461631e-13 13/35 156/9703
6 Atrial Fibrillation 3.624918e-13 13/35 160/9703
7 Auriculo-Ventricular Dissociation 4.078791e-03 2/35 4/9703
48 Heart Block 4.078791e-03 2/35 4/9703
185 Thyroid Agenesis 5.813634e-03 2/35 5/9703
41 Follicular cyst 2.597403e-02 1/35 1/9703
49 Cardiomegaly 2.597403e-02 3/35 82/9703
67 Malaria 2.597403e-02 2/35 20/9703******************************************
* *
* Welcome to WebGestaltR ! *
* *
******************************************
Loading the functional categories...
Loading the ID list...
Loading the reference list...
Performing the enrichment analysis...
description size overlap FDR database
1 Atrial Fibrillation 52 6 0.002615339 disease_GLAD4U
2 Arrhythmias, Cardiac 89 6 0.023624920 disease_GLAD4U
3 Isaacs Syndrome 55 5 0.023624920 disease_GLAD4U
userId
1 PRRX1;SCN10A;NKX2-5;CAV1;SYNPO2L;KCNJ5
2 PRRX1;SCN10A;NKX2-5;CAV1;SYNPO2L;KCNJ5
3 SCN10A;GNB4;NKX2-5;CCND2;KLF12PIP Manhattan Plot
library(tibble)
library(tidyverse)── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──✔ tidyr 1.1.4 ✔ dplyr 1.0.7
✔ readr 2.1.1 ✔ stringr 1.4.0
✔ purrr 0.3.4 ✔ forcats 0.5.1── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ tidyr::extract() masks disgenet2r::extract()
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()full.gene.pip.summary <- data.frame(gene_name = ctwas_gene_res$genename,
gene_pip = ctwas_gene_res$susie_pip,
gene_id = ctwas_gene_res$id,
chr = as.integer(ctwas_gene_res$chrom),
start = ctwas_gene_res$pos / 1e3,
is_highlight = F, stringsAsFactors = F) %>% as_tibble()
full.gene.pip.summary$is_highlight <- full.gene.pip.summary$gene_pip > 0.50
don <- full.gene.pip.summary %>%
# Compute chromosome size
group_by(chr) %>%
summarise(chr_len=max(start)) %>%
# Calculate cumulative position of each chromosome
mutate(tot=cumsum(chr_len)-chr_len) %>%
dplyr::select(-chr_len) %>%
# Add this info to the initial dataset
left_join(full.gene.pip.summary, ., by=c("chr"="chr")) %>%
# Add a cumulative position of each SNP
arrange(chr, start) %>%
mutate( BPcum=start+tot)
axisdf <- don %>% group_by(chr) %>% summarize(center=( max(BPcum) + min(BPcum) ) / 2 )
x_axis_labels <- axisdf$chr
x_axis_labels[seq(1,21,2)] <- ""
ggplot(don, aes(x=BPcum, y=gene_pip)) +
# Show all points
ggrastr::geom_point_rast(aes(color=as.factor(chr)), size=2) +
scale_color_manual(values = rep(c("grey", "skyblue"), 22 )) +
# custom X axis:
# scale_x_continuous(label = axisdf$chr,
# breaks= axisdf$center,
# guide = guide_axis(n.dodge = 2)) +
scale_x_continuous(label = x_axis_labels,
breaks = axisdf$center) +
scale_y_continuous(expand = c(0, 0), limits = c(0,1.25), breaks=(1:5)*0.2, minor_breaks=(1:10)*0.1) + # remove space between plot area and x axis
# Add highlighted points
ggrastr::geom_point_rast(data=subset(don, is_highlight==T), color="orange", size=2) +
# Add label using ggrepel to avoid overlapping
ggrepel::geom_label_repel(data=subset(don, is_highlight==T),
aes(label=gene_name),
size=4,
min.segment.length = 0,
label.size = NA,
fill = alpha(c("white"),0)) +
# Custom the theme:
theme_bw() +
theme(
text = element_text(size = 14),
legend.position="none",
panel.border = element_blank(),
panel.grid.major.x = element_blank(),
panel.grid.minor.x = element_blank()
) +
xlab("Chromosome") +
ylab("cTWAS PIP")Warning: ggrepel: 21 unlabeled data points (too many overlaps). Consider
increasing max.overlaps
| Version | Author | Date |
|---|---|---|
| 77de9fb | sq-96 | 2021-12-22 |
Locus Plots - 5_45
ctwas_gene_res_sortz <- ctwas_gene_res[order(-abs(ctwas_gene_res$z)),]
report_cols_region <- report_cols[!(report_cols %in% c("num_eqtl"))]locus_plot <- function(region_tag, rerun_ctwas = F, plot_eqtl = T, label="cTWAS"){
region_tag1 <- unlist(strsplit(region_tag, "_"))[1]
region_tag2 <- unlist(strsplit(region_tag, "_"))[2]
a <- ctwas_res[ctwas_res$region_tag==region_tag,]
regionlist <- readRDS(paste0(results_dir, "/", analysis_id, "_ctwas.regionlist.RDS"))
region <- regionlist[[as.numeric(region_tag1)]][[region_tag2]]
R_snp_info <- do.call(rbind, lapply(region$regRDS, function(x){data.table::fread(paste0(tools::file_path_sans_ext(x), ".Rvar"))}))
if (isTRUE(rerun_ctwas)){
ld_exprfs <- paste0(results_dir, "/", analysis_id, "_expr_chr", 1:22, ".expr.gz")
temp_reg <- data.frame("chr" = paste0("chr",region_tag1), "start" = region$start, "stop" = region$stop)
write.table(temp_reg,
#file= paste0(results_dir, "/", analysis_id, "_ctwas.temp.reg.txt") ,
file= "temp_reg.txt",
row.names=F, col.names=T, sep="\t", quote = F)
load(paste0(results_dir, "/", analysis_id, "_expr_z_snp.Rd"))
z_gene_temp <- z_gene[z_gene$id %in% a$id[a$type=="gene"],]
z_snp_temp <- z_snp[z_snp$id %in% R_snp_info$id,]
ctwas_rss(z_gene_temp, z_snp_temp, ld_exprfs, ld_pgenfs = NULL,
ld_R_dir = dirname(region$regRDS)[1],
ld_regions_custom = "temp_reg.txt", thin = 1,
outputdir = ".", outname = "temp", ncore = 1, ncore.rerun = 1, prob_single = 0,
group_prior = estimated_group_prior, group_prior_var = estimated_group_prior_var,
estimate_group_prior = F, estimate_group_prior_var = F)
a <- data.table::fread("temp.susieIrss.txt", header = T)
rownames(z_snp_temp) <- z_snp_temp$id
z_snp_temp <- z_snp_temp[a$id[a$type=="SNP"],]
z_gene_temp <- z_gene_temp[a$id[a$type=="gene"],]
a$z <- NA
a$z[a$type=="SNP"] <- z_snp_temp$z
a$z[a$type=="gene"] <- z_gene_temp$z
}
a$ifcausal <- 0
focus <- a$id[a$type=="gene"][which.max(abs(a$z[a$type=="gene"]))]
a$ifcausal <- as.numeric(a$id==focus)
a$PVALUE <- (-log(2) - pnorm(abs(a$z), lower.tail=F, log.p=T))/log(10)
R_gene <- readRDS(region$R_g_file)
R_snp_gene <- readRDS(region$R_sg_file)
R_snp <- as.matrix(Matrix::bdiag(lapply(region$regRDS, readRDS)))
rownames(R_gene) <- region$gid
colnames(R_gene) <- region$gid
rownames(R_snp_gene) <- R_snp_info$id
colnames(R_snp_gene) <- region$gid
rownames(R_snp) <- R_snp_info$id
colnames(R_snp) <- R_snp_info$id
a$r2max <- NA
a$r2max[a$type=="gene"] <- R_gene[focus,a$id[a$type=="gene"]]
a$r2max[a$type=="SNP"] <- R_snp_gene[a$id[a$type=="SNP"],focus]
r2cut <- 0.4
colorsall <- c("#7fc97f", "#beaed4", "#fdc086")
layout(matrix(1:2, ncol = 1), widths = 1, heights = c(1.5,1.5), respect = FALSE)
par(mar = c(0, 4.1, 4.1, 2.1))
plot(a$pos[a$type=="SNP"], a$PVALUE[a$type == "SNP"], pch = 19, xlab=paste0("Chromosome ", region_tag1, " Position"),frame.plot=FALSE, col = "white", ylim= c(-0.1,1.1), ylab = "cTWAS PIP", xaxt = 'n')
grid()
points(a$pos[a$type=="SNP"], a$susie_pip[a$type == "SNP"], pch = 21, xlab="Genomic position", bg = colorsall[1])
points(a$pos[a$type=="SNP" & a$r2max > r2cut], a$susie_pip[a$type == "SNP" & a$r2max >r2cut], pch = 21, bg = "purple")
points(a$pos[a$type=="SNP" & a$ifcausal == 1], a$susie_pip[a$type == "SNP" & a$ifcausal == 1], pch = 21, bg = "salmon")
points(a$pos[a$type=="gene"], a$susie_pip[a$type == "gene"], pch = 22, bg = colorsall[1], cex = 2)
points(a$pos[a$type=="gene" & a$r2max > r2cut], a$susie_pip[a$type == "gene" & a$r2max > r2cut], pch = 22, bg = "purple", cex = 2)
points(a$pos[a$type=="gene" & a$ifcausal == 1], a$susie_pip[a$type == "gene" & a$ifcausal == 1], pch = 22, bg = "salmon", cex = 2)
if (isTRUE(plot_eqtl)){
for (cgene in a[a$type=="gene" & a$ifcausal == 1, ]$id){
load(paste0(results_dir, "/",analysis_id, "_expr_chr", region_tag1, ".exprqc.Rd"))
eqtls <- rownames(wgtlist[[cgene]])
points(a[a$id %in% eqtls,]$pos, rep( -0.15, nrow(a[a$id %in% eqtls,])), pch = "|", col = "salmon", cex = 1.5)
}
}
legend(min(a$pos), y= 1.1 ,c("Gene", "SNP"), pch = c(22,21), title="Shape Legend", bty ='n', cex=0.6, title.adj = 0)
legend(min(a$pos), y= 0.7 ,c("Lead TWAS Gene", "R2 > 0.4", "R2 <= 0.4"), pch = 19, col = c("salmon", "purple", colorsall[1]), title="Color Legend", bty ='n', cex=0.6, title.adj = 0)
if (label=="cTWAS"){
text(a$pos[a$id==focus], a$susie_pip[a$id==focus], labels=ctwas_gene_res$genename[ctwas_gene_res$id==focus], pos=3, cex=0.6)
}
par(mar = c(4.1, 4.1, 0.5, 2.1))
plot(a$pos[a$type=="SNP"], a$PVALUE[a$type == "SNP"], pch = 21, xlab=paste0("Chromosome ", region_tag1, " Position"), frame.plot=FALSE, bg = colorsall[1], ylab = "TWAS -log10(p value)", panel.first = grid(), ylim =c(0, max(a$PVALUE)*1.2))
points(a$pos[a$type=="SNP" & a$r2max > r2cut], a$PVALUE[a$type == "SNP" & a$r2max > r2cut], pch = 21, bg = "purple")
points(a$pos[a$type=="SNP" & a$ifcausal == 1], a$PVALUE[a$type == "SNP" & a$ifcausal == 1], pch = 21, bg = "salmon")
points(a$pos[a$type=="gene"], a$PVALUE[a$type == "gene"], pch = 22, bg = colorsall[1], cex = 2)
points(a$pos[a$type=="gene" & a$r2max > r2cut], a$PVALUE[a$type == "gene" & a$r2max > r2cut], pch = 22, bg = "purple", cex = 2)
points(a$pos[a$type=="gene" & a$ifcausal == 1], a$PVALUE[a$type == "gene" & a$ifcausal == 1], pch = 22, bg = "salmon", cex = 2)
abline(h=-log10(alpha/nrow(ctwas_gene_res)), col ="red", lty = 2)
if (label=="TWAS"){
text(a$pos[a$id==focus], a$PVALUE[a$id==focus], labels=ctwas_gene_res$genename[ctwas_gene_res$id==focus], pos=3, cex=0.6)
}
}
locus_plot("10_49", label="TWAS")
| Version | Author | Date |
|---|---|---|
| 82c68fd | sq-96 | 2021-12-22 |
Genes with many eQTL
#distribution of number of eQTL for all imputed genes (after dropping ambiguous variants)
table(ctwas_gene_res$num_eqtl)
1 2 3 4 5 6
7076 3067 605 98 13 3 #all genes with 4+ eQTL
ctwas_gene_res[ctwas_gene_res$num_eqtl>3,] chrom id pos type region_tag1 region_tag2 cs_index
8829 1 ENSG00000169598.15 3856947 gene 1 3 0
7228 1 ENSG00000158825.5 20587958 gene 1 14 0
10362 1 ENSG00000183682.7 39491437 gene 1 24 0
11907 1 ENSG00000213366.12 109668049 gene 1 69 0
3657 1 ENSG00000120332.15 175068086 gene 1 86 0
13141 1 ENSG00000259865.1 247164089 gene 1 131 0
12651 2 ENSG00000242282.6 3525259 gene 2 2 0
13152 2 ENSG00000260077.1 10011416 gene 2 6 0
6061 2 ENSG00000144026.11 95125518 gene 2 57 0
3534 2 ENSG00000119147.9 106041343 gene 2 63 0
837 2 ENSG00000072163.19 127617901 gene 2 75 0
11022 2 ENSG00000196141.13 200299060 gene 2 118 0
5431 2 ENSG00000138363.14 215075224 gene 2 127 0
5865 2 ENSG00000142330.19 240586128 gene 2 143 0
11971 3 ENSG00000214021.15 9809721 gene 3 8 0
4655 4 ENSG00000132406.11 4247756 gene 4 5 0
12866 4 ENSG00000251129.1 31994228 gene 4 27 0
12748 4 ENSG00000246375.2 88284457 gene 4 59 0
7927 4 ENSG00000164124.10 158173634 gene 4 102 0
12863 4 ENSG00000250971.1 186813865 gene 4 120 0
12845 5 ENSG00000250490.1 6363492 gene 5 6 0
7951 5 ENSG00000164237.8 10308070 gene 5 9 0
4716 5 ENSG00000132840.9 79069112 gene 5 46 0
8063 5 ENSG00000164904.17 126544897 gene 5 77 0
12810 5 ENSG00000249306.5 174323562 gene 5 105 0
6251 6 ENSG00000145949.9 2680498 gene 6 3 0
6252 6 ENSG00000145979.17 13279232 gene 6 12 0
11690 6 ENSG00000204516.9 31494471 gene 6 27 0
11695 6 ENSG00000204531.17 31179417 gene 6 27 0
11944 6 ENSG00000213780.10 30904497 gene 6 27 0
11647 6 ENSG00000204301.6 32202656 gene 6 27 0
11633 6 ENSG00000204228.3 33204859 gene 6 27 0
9023 6 ENSG00000170915.8 52335661 gene 6 39 0
2358 7 ENSG00000106524.8 16599990 gene 7 16 0
6768 7 ENSG00000152926.14 64998115 gene 7 43 0
12433 7 ENSG00000234444.9 64307784 gene 7 43 0
11275 7 ENSG00000197558.11 149756134 gene 7 92 0
10781 7 ENSG00000187260.15 151404454 gene 7 94 0
2296 7 ENSG00000105983.20 156806908 gene 7 98 0
10215 8 ENSG00000182372.9 1292383 gene 8 3 0
12980 8 ENSG00000255394.4 11758186 gene 8 15 0
5316 8 ENSG00000137563.11 63029281 gene 8 48 0
11294 9 ENSG00000197646.7 5483361 gene 9 6 0
2416 9 ENSG00000107186.16 13106599 gene 9 12 0
5242 9 ENSG00000137074.18 33024919 gene 9 25 0
11958 9 ENSG00000213930.11 34637693 gene 9 27 0
6423 9 ENSG00000148357.16 130254358 gene 9 67 0
8177 9 ENSG00000165689.16 136406987 gene 9 73 0
5186 10 ENSG00000136738.14 17641079 gene 10 14 0
4791 10 ENSG00000133661.15 79979389 gene 10 51 0
5404 10 ENSG00000138119.16 93323282 gene 10 59 0
8952 10 ENSG00000170430.9 129467281 gene 10 81 0
9724 11 ENSG00000177042.14 688091 gene 11 1 0
3726 11 ENSG00000121236.20 5590335 gene 11 4 0
2696 11 ENSG00000110328.5 11410495 gene 11 9 0
764 12 ENSG00000069493.14 9664279 gene 12 9 0
40 12 ENSG00000004700.15 21501331 gene 12 16 0
11169 12 ENSG00000196935.8 63844067 gene 12 39 0
2808 12 ENSG00000111581.9 68686129 gene 12 42 0
10495 12 ENSG00000184967.6 132083372 gene 12 81 0
13025 12 ENSG00000256576.2 132189571 gene 12 81 0
2018 15 ENSG00000103966.10 41923122 gene 15 15 0
10624 15 ENSG00000185880.12 44727667 gene 15 17 0
1214 16 ENSG00000086504.15 345757 gene 16 1 0
1938 16 ENSG00000103148.15 137677 gene 16 1 0
1980 16 ENSG00000103381.11 12803538 gene 16 13 0
10341 16 ENSG00000183549.10 20409610 gene 16 19 0
6544 16 ENSG00000149922.10 30091070 gene 16 24 0
5732 16 ENSG00000140955.10 84184959 gene 16 48 0
5743 16 ENSG00000141013.16 90018731 gene 16 54 0
4344 17 ENSG00000129204.16 5116034 gene 17 5 0
10120 17 ENSG00000181291.6 34581017 gene 17 21 0
13108 17 ENSG00000259207.7 47221042 gene 17 27 0
13351 17 ENSG00000266714.6 75588005 gene 17 42 0
13289 17 ENSG00000262877.4 81388422 gene 17 46 0
5796 17 ENSG00000141526.16 82228706 gene 17 47 0
944 18 ENSG00000075643.5 36187019 gene 18 19 0
708 18 ENSG00000066926.10 57571588 gene 18 31 0
651 19 ENSG00000065268.10 982793 gene 19 2 0
1526 19 ENSG00000099817.11 1032229 gene 19 2 0
976 19 ENSG00000076826.9 7594599 gene 19 8 0
277 19 ENSG00000021488.12 32869435 gene 19 23 0
13395 19 ENSG00000267475.1 32691414 gene 19 23 0
11312 19 ENSG00000197782.14 40070352 gene 19 27 0
4459 19 ENSG00000130529.15 49114644 gene 19 35 0
8552 19 ENSG00000167766.18 52615913 gene 19 36 0
9027 19 ENSG00000170949.17 53087254 gene 19 36 0
11285 20 ENSG00000197586.12 25195384 gene 20 18 0
10660 20 ENSG00000186191.7 33077982 gene 20 20 0
11015 20 ENSG00000196090.12 43187301 gene 20 27 0
12140 20 ENSG00000223891.5 44210535 gene 20 28 0
7376 21 ENSG00000160305.17 46438335 gene 21 24 0
1440 22 ENSG00000093072.15 17197305 gene 22 2 0
1546 22 ENSG00000099957.16 21003838 gene 22 4 0
10551 22 ENSG00000185339.8 30596847 gene 22 10 0
10343 22 ENSG00000183569.17 42552989 gene 22 18 0
10997 22 ENSG00000189306.10 42490805 gene 22 18 0
10747 22 ENSG00000186976.14 43808777 gene 22 19 0
99 2 ENSG00000006607.13 241354886 gene 2 144 0
8665 2 ENSG00000168397.16 241636141 gene 2 144 0
10082 2 ENSG00000180902.17 241732627 gene 2 144 0
8101 6 ENSG00000165097.13 18153891 gene 6 14 0
8296 10 ENSG00000166295.8 72211079 gene 10 49 0
9454 11 ENSG00000174370.9 128809746 gene 11 80 0
6594 12 ENSG00000150977.10 123430286 gene 12 75 0
7162 12 ENSG00000158023.9 121918251 gene 12 75 0
3866 13 ENSG00000123179.13 49691102 gene 13 21 0
7111 14 ENSG00000157379.13 24291558 gene 14 3 0
5639 14 ENSG00000140043.11 73844725 gene 14 34 0
163 14 ENSG00000009830.11 77311261 gene 14 36 0
10088 15 ENSG00000180953.11 79922366 gene 15 37 0
1334 19 ENSG00000089847.12 4183158 gene 19 4 0
7353 21 ENSG00000160200.17 42941853 gene 21 22 0
7358 21 ENSG00000160213.6 43767678 gene 21 22 0
susie_pip mu2 region_tag PVE genename
8829 4.160727e-02 4.562732 1_3 1.841639e-07 DFFB
7228 1.659371e-01 16.299471 1_14 2.623781e-06 CDA
10362 3.536819e-02 4.661882 1_24 1.599501e-07 BMP8A
11907 1.142417e-02 9.805167 1_69 1.086651e-07 GSTM2
3657 6.710729e-02 12.989658 1_86 8.456251e-07 TNN
13141 3.085086e-02 9.011259 1_131 2.696890e-07 RP11-488L18.10
12651 1.896278e-02 7.677279 2_2 1.412277e-07 AC108488.4
13152 4.407680e-02 6.839611 2_6 2.924502e-07 RP11-254F7.2
6061 8.657371e-02 8.343336 2_57 7.007066e-07 ZNF514
3534 1.411348e-01 14.521299 2_63 1.988153e-06 C2orf40
837 6.084210e-02 5.313185 2_75 3.135953e-07 LIMS2
11022 3.003384e-02 95.300331 2_118 2.776615e-06 SPATS2L
5431 1.360629e-01 16.535574 2_127 2.182576e-06 ATIC
5865 3.482796e-02 7.164635 2_143 2.420653e-07 CAPN10
11971 1.188509e-01 11.239305 3_8 1.295843e-06 TTLL3
4655 4.169085e-02 4.603728 4_5 1.861919e-07 TMEM128
12866 6.874848e-02 7.086021 4_27 4.725806e-07 RP11-734I18.1
12748 5.990045e-02 4.698849 4_59 2.730436e-07 RP11-10L7.1
7927 6.030136e-02 4.990715 4_102 2.919445e-07 TMEM144
12863 6.017045e-02 5.706989 4_120 3.331200e-07 RP11-696F12.1
12845 1.198806e-01 13.958888 5_6 1.623343e-06 LINC02145
7951 5.059846e-02 5.388011 5_9 2.644698e-07 CMBL
4716 4.239101e-02 5.025132 5_46 2.066482e-07 BHMT2
8063 6.077566e-02 5.731451 5_77 3.379128e-07 ALDH7A1
12810 6.464465e-02 7.490826 5_105 4.697564e-07 LINC01411
6251 4.335437e-02 4.532337 6_3 1.906187e-07 MYLK4
6252 5.601621e-02 5.601698 6_12 3.043994e-07 TBC1D7
11690 6.007886e-02 27.777389 6_27 1.618913e-06 MICB
11695 8.047401e-03 9.415482 6_27 7.350360e-08 POU5F1
11944 1.296462e-02 13.811198 6_27 1.737007e-07 GTF2H4
11647 4.710762e-03 4.543353 6_27 2.076243e-08 NOTCH4
11633 1.396179e-02 14.381660 6_27 1.947873e-07 HSD17B8
9023 2.867181e-02 4.543987 6_39 1.263871e-07 PAQR8
2358 4.702895e-02 5.445215 7_16 2.484224e-07 ANKMY2
6768 5.023305e-02 6.151104 7_43 2.997458e-07 ZNF117
12433 8.419199e-02 11.240276 7_43 9.180328e-07 ZNF736
11275 5.200805e-02 8.969464 7_92 4.525301e-07 SSPO
10781 3.396058e-02 8.645892 7_94 2.848363e-07 WDR86
2296 6.147956e-02 8.163856 7_98 4.868963e-07 LMBR1
10215 4.632397e-02 4.535147 8_3 2.038016e-07 CLN8
12980 2.476107e-02 8.773421 8_15 2.107408e-07 C8orf49
5316 6.841080e-02 7.031033 8_48 4.666102e-07 GGH
11294 5.306128e-02 4.855033 9_6 2.499081e-07 PDCD1LG2
2416 3.206349e-02 4.883384 9_12 1.518945e-07 MPDZ
5242 5.319653e-02 4.542096 9_25 2.343959e-07 APTX
11958 2.328524e-02 5.078618 9_27 1.147193e-07 GALT
6423 2.466585e-01 20.836963 9_67 4.985871e-06 HMCN2
8177 8.784738e-02 13.072608 9_73 1.114042e-06 SDCCAG3
5186 5.387607e-02 10.479376 10_14 5.476987e-07 STAM
4791 6.836780e-02 8.431123 10_51 5.591746e-07 SFTPD
5404 1.445696e-01 15.503954 10_59 2.174352e-06 MYOF
8952 6.789768e-02 9.622812 10_81 6.338221e-07 MGMT
9724 3.960285e-02 4.552221 11_1 1.748881e-07 TMEM80
3726 3.951884e-02 4.530694 11_4 1.736918e-07 TRIM6
2696 4.172490e-02 4.834936 11_9 1.957025e-07 GALNT18
764 6.774877e-02 6.210837 12_9 4.081896e-07 CLEC2D
40 4.723176e-02 11.482884 12_16 5.261329e-07 RECQL
11169 5.232829e-02 6.129059 12_39 3.111292e-07 SRGAP1
2808 1.228938e-01 13.104567 12_42 1.562295e-06 NUP107
10495 1.028720e-01 12.642046 12_81 1.261609e-06 NOC4L
13025 4.396147e-02 4.711916 12_81 2.009464e-07 RP13-977J11.2
2018 1.444954e-01 17.024716 15_15 2.386407e-06 EHD4
10624 6.678218e-02 5.884794 15_17 3.812433e-07 TRIM69
1214 4.551849e-02 8.109793 16_1 3.581031e-07 MRPL28
1938 1.796355e-01 21.054359 16_1 3.668973e-06 NPRL3
1980 5.161597e-02 6.782632 16_13 3.396196e-07 CPPED1
10341 9.125280e-02 9.051325 16_19 8.012514e-07 ACSM5
6544 4.402615e-02 5.613910 16_24 2.397654e-07 TBX6
5732 6.711738e-02 7.912542 16_48 5.151830e-07 ADAD2
5743 5.375185e-02 6.618971 16_54 3.451392e-07 GAS8
4344 4.956765e-02 5.281360 17_5 2.539536e-07 USP6
10120 2.682028e-02 4.984821 17_21 1.296950e-07 TMEM132E
13108 1.319177e-05 10.880709 17_27 1.392421e-10 ITGB3
13351 3.910169e-02 4.592795 17_42 1.742140e-07 MYO15B
13289 4.163151e-02 5.107075 17_46 2.062552e-07 RP11-1055B8.4
5796 5.954408e-02 4.816236 17_47 2.781998e-07 SLC16A3
944 4.326557e-02 4.562763 18_19 1.915053e-07 MOCOS
708 5.743433e-02 7.227610 18_31 4.026954e-07 FECH
651 4.324470e-02 7.698742 19_2 3.229707e-07 WDR18
1526 2.684089e-01 25.177110 19_2 6.555612e-06 POLR2E
976 1.960340e-02 4.535871 19_8 8.625865e-08 CAMSAP3
277 1.329474e-01 17.026823 19_23 2.195957e-06 SLC7A9
13395 3.501099e-02 4.556265 19_23 1.547476e-07 CTD-2538C1.2
11312 3.907937e-02 4.680517 19_27 1.774401e-07 ZNF780A
4459 1.348543e-02 5.099611 19_35 6.671329e-08 TRPM4
8552 3.337832e-02 4.565379 19_36 1.478263e-07 ZNF83
9027 3.869678e-02 5.924854 19_36 2.224144e-07 ZNF160
11285 5.941315e-02 8.129045 20_18 4.685247e-07 ENTPD6
10660 5.525131e-02 6.947800 20_20 3.723920e-07 BPIFB4
11015 3.867975e-02 4.572127 20_27 1.715586e-07 PTPRT
12140 2.660360e-02 6.086123 20_28 1.570694e-07 OSER1-AS1
7376 3.970524e-02 6.138704 21_24 2.364476e-07 DIP2A
1440 4.014696e-02 4.625767 22_2 1.801552e-07 CECR1
1546 1.909623e-02 7.224043 22_4 1.338253e-07 P2RX6
10551 1.456716e-01 16.721527 22_10 2.362986e-06 TCN2
10343 4.385792e-02 4.696390 22_18 1.998125e-07 SERHL2
10997 4.472207e-02 4.876429 22_18 2.115603e-07 RRP7A
10747 5.866690e-02 6.643992 22_19 3.781226e-07 EFCAB6
99 3.773174e-02 4.920054 2_144 1.800889e-07 FARP2
8665 3.819760e-02 5.543522 2_144 2.054150e-07 ATG4B
10082 3.665308e-02 4.592015 2_144 1.632767e-07 D2HGDH
8101 3.367752e-02 47.184002 6_14 1.541506e-06 KDM1B
8296 7.630441e-05 11.415357 10_49 8.449861e-10 ANAPC16
9454 2.571903e-02 19.573814 11_80 4.883605e-07 C11orf45
6594 9.310735e-04 4.788344 12_75 4.324937e-09 RILPL2
7162 1.577479e-03 9.393574 12_75 1.437490e-08 WDR66
3866 5.454124e-02 10.259270 13_21 5.428150e-07 EBPL
7111 1.370814e-02 5.168687 14_3 6.873363e-08 DHRS1
5639 1.457521e-02 5.222525 14_34 7.384237e-08 PTGR2
163 2.537794e-02 4.919146 14_36 1.211034e-07 POMT2
10088 2.986045e-02 7.646866 15_37 2.215084e-07 ST20
1334 3.238781e-02 7.550345 19_4 2.372241e-07 ANKRD24
7353 1.054189e-02 5.401357 21_22 5.523721e-08 CBS
7358 1.076657e-02 5.588965 21_22 5.837398e-08 CSTB
gene_type z num_eqtl
8829 protein_coding -0.05606883 4
7228 protein_coding -2.11159404 5
10362 protein_coding -0.04886818 4
11907 protein_coding -1.13618389 4
3657 protein_coding -1.56712336 4
13141 lincRNA -1.15241658 4
12651 lincRNA 0.88509724 4
13152 lincRNA -0.65540955 5
6061 protein_coding 1.18030677 5
3534 protein_coding 1.75470100 5
837 protein_coding 0.67333098 4
11022 protein_coding 8.11412691 4
5431 protein_coding -1.95877307 4
5865 protein_coding -1.08230144 4
11971 protein_coding 1.41366703 4
4655 protein_coding 0.35456062 4
12866 lincRNA 0.99806910 4
12748 lincRNA 0.16815663 4
7927 protein_coding -0.36615299 5
12863 lincRNA 0.52938089 4
12845 lincRNA -1.84178478 4
7951 protein_coding 0.57035749 4
4716 protein_coding 0.55240521 4
8063 protein_coding 0.67205955 4
12810 lincRNA -0.87537902 4
6251 protein_coding -0.01538851 4
6252 protein_coding 0.58335172 4
11690 protein_coding 2.73153856 4
11695 protein_coding 1.19707710 4
11944 protein_coding 1.54571781 4
11647 protein_coding -0.04864491 4
11633 protein_coding -1.70351244 4
9023 protein_coding 0.05092472 4
2358 protein_coding 0.59449010 4
6768 protein_coding 0.29410631 4
12433 protein_coding -1.42962198 4
11275 protein_coding 1.17327949 4
10781 protein_coding -1.07155972 4
2296 protein_coding 1.12846744 4
10215 protein_coding 0.14155900 5
12980 lincRNA 2.00933006 4
5316 protein_coding -1.03772803 5
11294 protein_coding -0.13876175 4
2416 protein_coding -0.31465621 4
5242 protein_coding -0.10240511 4
11958 protein_coding 0.54117071 4
6423 protein_coding 2.58282602 6
8177 protein_coding 1.62552220 4
5186 protein_coding -1.21519734 5
4791 protein_coding -1.06876491 4
5404 protein_coding 1.89158295 4
8952 protein_coding -1.23104737 4
9724 protein_coding -0.21243537 4
3726 protein_coding 0.03791292 4
2696 protein_coding -0.25931605 4
764 protein_coding -0.74959674 4
40 protein_coding 1.39889693 4
11169 protein_coding -0.64611366 4
2808 protein_coding 1.57728136 4
10495 protein_coding -1.59211541 6
13025 lincRNA -0.17676741 4
2018 protein_coding -2.18379904 4
10624 protein_coding 0.68260923 4
1214 protein_coding 1.46503375 4
1938 protein_coding -2.38794751 5
1980 protein_coding -0.78395748 4
10341 protein_coding -1.23454259 4
6544 protein_coding 0.23682640 4
5732 protein_coding -1.14113140 4
5743 protein_coding -0.82726355 4
4344 protein_coding -0.28654155 4
10120 protein_coding -0.33913678 4
13108 protein_coding -1.57910576 4
13351 protein_coding 0.07613186 5
13289 lincRNA -0.44622658 4
5796 protein_coding -0.36078965 4
944 protein_coding -0.13469930 4
708 protein_coding 0.79795344 4
651 protein_coding 1.06428772 4
1526 protein_coding -2.83809922 4
976 protein_coding -0.01811536 4
277 protein_coding -2.02040773 4
13395 lincRNA 0.08674841 4
11312 protein_coding 0.34051207 4
4459 protein_coding -0.38500405 4
8552 protein_coding -0.10819463 5
9027 protein_coding -0.58121153 4
11285 protein_coding 1.11446112 4
10660 protein_coding 0.88525055 4
11015 protein_coding 0.06492351 4
12140 lincRNA -0.42101203 5
7376 protein_coding 0.60857479 4
1440 protein_coding -0.16372486 4
1546 protein_coding -0.82364067 4
10551 protein_coding 2.42724084 4
10343 protein_coding 0.22201409 4
10997 protein_coding -0.17607349 6
10747 protein_coding -0.82194388 4
99 protein_coding 0.51764732 4
8665 protein_coding 1.00584632 4
10082 protein_coding 0.02223367 4
8101 protein_coding -7.80377324 4
8296 protein_coding 1.50960211 4
9454 protein_coding -4.26886900 4
6594 protein_coding -0.31787410 4
7162 protein_coding -1.03591335 4
3866 protein_coding -0.71354410 4
7111 protein_coding 0.28913514 4
5639 protein_coding -0.74125560 4
163 protein_coding -0.09069729 5
10088 protein_coding 1.24686579 4
1334 protein_coding -0.95784361 4
7353 protein_coding -0.47878923 4
7358 protein_coding -0.62864698 4#distribution of number of eQTL for genes with PIP>0.8
table(ctwas_gene_res$num_eqtl[ctwas_gene_res$susie_pip>0.8])/sum(ctwas_gene_res$susie_pip>0.8)
1 2 3
0.625 0.250 0.125 #genes with 2+ eQTL and PIP>0.8
ctwas_gene_res[ctwas_gene_res$num_eqtl>1 & ctwas_gene_res$susie_pip>0.8,] chrom id pos type region_tag1 region_tag2 cs_index
3275 1 ENSG00000116132.11 170662622 gene 1 84 4
6114 2 ENSG00000144589.21 219597759 gene 2 130 0
712 2 ENSG00000067066.16 230415508 gene 2 135 0
9691 2 ENSG00000176720.4 241555444 gene 2 144 0
2293 7 ENSG00000105974.11 116519968 gene 7 70 1
9012 8 ENSG00000170873.18 124576649 gene 8 82 1
8992 9 ENSG00000170681.6 100578027 gene 9 50 0
9257 10 ENSG00000172650.13 73650294 gene 10 49 5
2444 10 ENSG00000107651.12 119893621 gene 10 74 1
11818 14 ENSG00000205683.11 72894134 gene 14 34 1
2138 19 ENSG00000104964.14 3056215 gene 19 4 0
6914 21 ENSG00000154721.14 25638800 gene 21 9 1
susie_pip mu2 region_tag PVE genename gene_type
3275 0.9998808 119.34663 1_84 1.157627e-04 PRRX1 protein_coding
6114 0.8406381 18.96825 2_130 1.546845e-05 STK11IP protein_coding
712 0.8658723 18.73540 2_135 1.573719e-05 SP100 protein_coding
9691 0.8255975 19.18259 2_144 1.536335e-05 BOK protein_coding
2293 1.0000000 622.84815 7_70 6.042165e-04 CAV1 protein_coding
9012 0.8594108 20.87861 8_82 1.740655e-05 MTSS1 protein_coding
8992 0.8478869 23.34736 9_50 1.920376e-05 MURC protein_coding
9257 0.9779202 48.91516 10_49 4.640420e-05 AGAP5 protein_coding
2444 0.9517809 22.38862 10_74 2.067163e-05 SEC23IP protein_coding
11818 0.9574312 33.35512 14_34 3.097993e-05 DPF3 protein_coding
2138 0.9403219 20.20438 19_4 1.843030e-05 AES protein_coding
6914 0.9640543 22.19279 21_9 2.075505e-05 JAM2 protein_coding
z num_eqtl
3275 14.667578 2
6114 -3.868022 2
712 -3.671335 2
9691 3.910125 3
2293 15.567870 3
9012 4.402634 2
8992 4.911964 2
9257 11.518590 2
2444 -4.565228 2
11818 6.264960 3
2138 4.182804 3
6914 4.563232 2
sessionInfo()R version 3.6.1 (2019-07-05)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Scientific Linux 7.4 (Nitrogen)
Matrix products: default
BLAS/LAPACK: /software/openblas-0.2.19-el7-x86_64/lib/libopenblas_haswellp-r0.2.19.so
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] forcats_0.5.1 stringr_1.4.0 dplyr_1.0.7 purrr_0.3.4
[5] readr_2.1.1 tidyr_1.1.4 tidyverse_1.3.1 tibble_3.1.6
[9] WebGestaltR_0.4.4 disgenet2r_0.99.2 enrichR_3.0 cowplot_1.0.0
[13] ggplot2_3.3.5 workflowr_1.6.2
loaded via a namespace (and not attached):
[1] fs_1.5.2 lubridate_1.8.0 bit64_4.0.5 doParallel_1.0.16
[5] httr_1.4.2 rprojroot_2.0.2 tools_3.6.1 backports_1.4.1
[9] doRNG_1.8.2 utf8_1.2.2 R6_2.5.1 vipor_0.4.5
[13] DBI_1.1.1 colorspace_2.0-2 withr_2.4.3 ggrastr_1.0.1
[17] tidyselect_1.1.1 bit_4.0.4 curl_4.3.2 compiler_3.6.1
[21] git2r_0.26.1 cli_3.1.0 rvest_1.0.2 Cairo_1.5-12.2
[25] xml2_1.3.3 labeling_0.4.2 scales_1.1.1 apcluster_1.4.8
[29] digest_0.6.29 rmarkdown_2.11 svglite_1.2.2 pkgconfig_2.0.3
[33] htmltools_0.5.2 dbplyr_2.1.1 fastmap_1.1.0 highr_0.9
[37] rlang_0.4.12 readxl_1.3.1 rstudioapi_0.13 RSQLite_2.2.8
[41] jquerylib_0.1.4 farver_2.1.0 generics_0.1.1 jsonlite_1.7.2
[45] vroom_1.5.7 magrittr_2.0.1 Matrix_1.2-18 ggbeeswarm_0.6.0
[49] Rcpp_1.0.7 munsell_0.5.0 fansi_0.5.0 gdtools_0.1.9
[53] lifecycle_1.0.1 stringi_1.7.6 whisker_0.3-2 yaml_2.2.1
[57] plyr_1.8.6 grid_3.6.1 blob_1.2.2 ggrepel_0.9.1
[61] parallel_3.6.1 promises_1.0.1 crayon_1.4.2 lattice_0.20-38
[65] haven_2.4.3 hms_1.1.1 knitr_1.36 pillar_1.6.4
[69] igraph_1.2.10 rjson_0.2.20 rngtools_1.5.2 reshape2_1.4.4
[73] codetools_0.2-16 reprex_2.0.1 glue_1.5.1 evaluate_0.14
[77] data.table_1.14.2 modelr_0.1.8 vctrs_0.3.8 tzdb_0.2.0
[81] httpuv_1.5.1 foreach_1.5.1 cellranger_1.1.0 gtable_0.3.0
[85] assertthat_0.2.1 cachem_1.0.6 xfun_0.29 broom_0.7.10
[89] later_0.8.0 iterators_1.0.13 beeswarm_0.2.3 memoise_2.0.1
[93] ellipsis_0.3.2