Last updated: 2021-09-09
Checks: 6 1
Knit directory: ctwas_applied/
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.
The R Markdown file has unstaged changes. To know which version of the R Markdown file created these results, you’ll want to first commit it to the Git repo. If you’re still working on the analysis, you can ignore this warning. When you’re finished, you can run wflow_publish
to commit the R Markdown file and build the HTML.
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(20210726)
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.
Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.
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 59e5f4d. 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:
Unstaged changes:
Modified: analysis/ukb-d-30500_irnt_Liver.Rmd
Modified: analysis/ukb-d-30500_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30600_irnt_Liver.Rmd
Modified: analysis/ukb-d-30600_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30610_irnt_Liver.Rmd
Modified: analysis/ukb-d-30610_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30620_irnt_Liver.Rmd
Modified: analysis/ukb-d-30620_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30630_irnt_Liver.Rmd
Modified: analysis/ukb-d-30630_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30640_irnt_Liver.Rmd
Modified: analysis/ukb-d-30640_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30650_irnt_Liver.Rmd
Modified: analysis/ukb-d-30650_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30660_irnt_Liver.Rmd
Modified: analysis/ukb-d-30660_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30670_irnt_Liver.Rmd
Modified: analysis/ukb-d-30670_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30680_irnt_Liver.Rmd
Modified: analysis/ukb-d-30690_irnt_Liver.Rmd
Modified: analysis/ukb-d-30690_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30700_irnt_Liver.Rmd
Modified: analysis/ukb-d-30700_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30710_irnt_Liver.Rmd
Modified: analysis/ukb-d-30710_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30720_irnt_Liver.Rmd
Modified: analysis/ukb-d-30720_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30730_irnt_Liver.Rmd
Modified: analysis/ukb-d-30740_irnt_Liver.Rmd
Modified: analysis/ukb-d-30740_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30750_irnt_Liver.Rmd
Modified: analysis/ukb-d-30750_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30760_irnt_Liver.Rmd
Modified: analysis/ukb-d-30760_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30770_irnt_Liver.Rmd
Modified: analysis/ukb-d-30770_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30780_irnt_Liver.Rmd
Modified: analysis/ukb-d-30780_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30790_irnt_Liver.Rmd
Modified: analysis/ukb-d-30800_irnt_Liver.Rmd
Modified: analysis/ukb-d-30800_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30810_irnt_Liver.Rmd
Modified: analysis/ukb-d-30820_irnt_Liver.Rmd
Modified: analysis/ukb-d-30820_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30830_irnt_Liver.Rmd
Modified: analysis/ukb-d-30830_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30840_irnt_Liver.Rmd
Modified: analysis/ukb-d-30840_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30850_irnt_Liver.Rmd
Modified: analysis/ukb-d-30850_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30860_irnt_Liver.Rmd
Modified: analysis/ukb-d-30860_irnt_Whole_Blood.Rmd
Modified: analysis/ukb-d-30870_irnt_Liver.Rmd
Modified: analysis/ukb-d-30880_irnt_Liver.Rmd
Modified: analysis/ukb-d-30890_irnt_Liver.Rmd
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/ukb-d-30860_irnt_Whole_Blood.Rmd
) and HTML (docs/ukb-d-30860_irnt_Whole_Blood.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 | cbf7408 | wesleycrouse | 2021-09-08 | adding enrichment to reports |
html | cbf7408 | wesleycrouse | 2021-09-08 | adding enrichment to reports |
Rmd | 4970e3e | wesleycrouse | 2021-09-08 | updating reports |
html | 4970e3e | wesleycrouse | 2021-09-08 | updating reports |
Rmd | dfd2b5f | wesleycrouse | 2021-09-07 | regenerating reports |
html | dfd2b5f | wesleycrouse | 2021-09-07 | regenerating reports |
Rmd | 61b53b3 | wesleycrouse | 2021-09-06 | updated PVE calculation |
html | 61b53b3 | wesleycrouse | 2021-09-06 | updated PVE calculation |
Rmd | 837dd01 | wesleycrouse | 2021-09-01 | adding additional fixedsigma report |
Rmd | d0a5417 | wesleycrouse | 2021-08-30 | adding new reports to the index |
Rmd | 0922de7 | wesleycrouse | 2021-08-18 | updating all reports with locus plots |
html | 0922de7 | wesleycrouse | 2021-08-18 | updating all reports with locus plots |
html | 1c62980 | wesleycrouse | 2021-08-11 | Updating reports |
Rmd | 5981e80 | wesleycrouse | 2021-08-11 | Adding more reports |
html | 5981e80 | wesleycrouse | 2021-08-11 | Adding more reports |
Rmd | da9f015 | wesleycrouse | 2021-08-07 | adding more reports |
html | da9f015 | wesleycrouse | 2021-08-07 | adding more reports |
These are the results of a ctwas
analysis of the UK Biobank trait Total protein (quantile)
using Whole_Blood
gene weights.
The GWAS was conducted by the Neale Lab, and the biomarker traits we analyzed are discussed here. Summary statistics were obtained from IEU OpenGWAS using GWAS ID: ukb-d-30860_irnt
. Results were obtained from from IEU rather than Neale Lab because they are in a standardard format (GWAS VCF). Note that 3 of the 34 biomarker traits were not available from IEU and were excluded from analysis.
The weights are mashr GTEx v8 models on Whole_Blood
eQTL obtained from PredictDB. We performed a full harmonization of the variants, including recovering strand ambiguous variants. This procedure is discussed in a separate document. (TO-DO: add report that describes harmonization)
LD matrices were computed from a 10% subset of Neale lab subjects. Subjects were matched using the plate and well information from genotyping. We included only biallelic variants with MAF>0.01 in the original Neale Lab GWAS. (TO-DO: add more details [number of subjects, variants, etc])
TO-DO: add enhanced QC reporting (total number of weights, why each variant was missing for all genes)
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] 11095
#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
1129 747 624 400 479 621 560 383 404 430 682 652 192 362 331
16 17 18 19 20 21 22
551 725 159 911 313 130 310
#proportion of imputed weights without missing variants
mean(qclist_all$nmiss==0)
[1] 0.762776
#load ctwas results
ctwas_res <- data.table::fread(paste0(results_dir, "/", analysis_id, "_ctwas.susieIrss.txt"))
#make unique identifier for regions
ctwas_res$region_tag <- paste(ctwas_res$region_tag1, ctwas_res$region_tag2, sep="_")
#compute PVE for each gene/SNP
ctwas_res$PVE = ctwas_res$susie_pip*ctwas_res$mu2/sample_size #check PVE calculation
#separate gene and SNP results
ctwas_gene_res <- ctwas_res[ctwas_res$type == "gene", ]
ctwas_gene_res <- data.frame(ctwas_gene_res)
ctwas_snp_res <- ctwas_res[ctwas_res$type == "SNP", ]
ctwas_snp_res <- data.frame(ctwas_snp_res)
#add gene information to results
sqlite <- RSQLite::dbDriver("SQLite")
db = RSQLite::dbConnect(sqlite, paste0("/project2/compbio/predictdb/mashr_models/mashr_", weight, ".db"))
query <- function(...) RSQLite::dbGetQuery(db, ...)
gene_info <- query("select gene, genename from extra")
gene_info <- query("select gene, genename, gene_type from extra")
RSQLite::dbDisconnect(db)
ctwas_gene_res <- cbind(ctwas_gene_res, gene_info[sapply(ctwas_gene_res$id, match, gene_info$gene), c("genename", "gene_type")])
#add z score to results
load(paste0(results_dir, "/", analysis_id, "_expr_z_gene.Rd"))
ctwas_gene_res$z <- z_gene[ctwas_gene_res$id,]$z
#load(paste0(results_dir, "/", analysis_id, "_expr_z_snp.Rd")) #for new version, stored after harmonization
z_snp <- readRDS(paste0(results_dir, "/", trait_id, ".RDS")) #for old version, unharmonized
z_snp <- z_snp[z_snp$id %in% ctwas_snp_res$id,] #subset snps to those included in analysis, note some are duplicated, need to match which allele was used
ctwas_snp_res$z <- z_snp$z[match(ctwas_snp_res$id, z_snp$id)] #for duplicated snps, this arbitrarily uses the first allele
ctwas_snp_res$z_flag <- as.numeric(ctwas_snp_res$id %in% z_snp$id[duplicated(z_snp$id)]) #mark the unclear z scores, flag=1
#formatting and rounding for tables
ctwas_gene_res$z <- round(ctwas_gene_res$z,2)
ctwas_snp_res$z <- round(ctwas_snp_res$z,2)
ctwas_gene_res$susie_pip <- round(ctwas_gene_res$susie_pip,3)
ctwas_snp_res$susie_pip <- round(ctwas_snp_res$susie_pip,3)
ctwas_gene_res$mu2 <- round(ctwas_gene_res$mu2,2)
ctwas_snp_res$mu2 <- round(ctwas_snp_res$mu2,2)
ctwas_gene_res$PVE <- signif(ctwas_gene_res$PVE, 2)
ctwas_snp_res$PVE <- signif(ctwas_snp_res$PVE, 2)
#merge gene and snp results with added information
ctwas_gene_res$z_flag=NA
ctwas_snp_res$genename=NA
ctwas_snp_res$gene_type=NA
ctwas_res <- rbind(ctwas_gene_res,
ctwas_snp_res[,colnames(ctwas_gene_res)])
#store columns to report
report_cols <- colnames(ctwas_gene_res)[!(colnames(ctwas_gene_res) %in% c("type", "region_tag1", "region_tag2", "cs_index", "gene_type", "z_flag", "id", "chrom", "pos"))]
first_cols <- c("genename", "region_tag")
report_cols <- c(first_cols, report_cols[!(report_cols %in% first_cols)])
report_cols_snps <- c("id", report_cols[-1])
#get number of SNPs from s1 results; adjust for thin
ctwas_res_s1 <- data.table::fread(paste0(results_dir, "/", analysis_id, "_ctwas.s1.susieIrss.txt"))
n_snps <- sum(ctwas_res_s1$type=="SNP")/thin
rm(ctwas_res_s1)
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)
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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.0137454575 0.0001980013
#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
25.87054 13.58456
#report sample size
print(sample_size)
[1] 314921
#report group size
group_size <- c(nrow(ctwas_gene_res), n_snps)
print(group_size)
[1] 11095 8697330
#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.01252825 0.07428446
#compare sum(PIP*mu2/sample_size) with above PVE calculation
c(sum(ctwas_gene_res$PVE),sum(ctwas_snp_res$PVE))
[1] 0.7911534 2.8245965
#distribution of PIPs
hist(ctwas_gene_res$susie_pip, xlim=c(0,1), main="Distribution of Gene PIPs")
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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
1980 FCGRT 19_34 1.000 238599.30 7.6e-01 17.84
6089 FADS1 11_34 0.997 187.54 5.9e-04 13.66
1541 ASCC2 22_10 0.997 662.80 2.1e-03 3.15
6824 TNFRSF13C 22_17 0.992 60.24 1.9e-04 7.57
5942 CHMP7 8_24 0.988 33.33 1.0e-04 -5.46
5631 ATP8B2 1_75 0.986 74.88 2.3e-04 9.37
7460 ANKRD55 5_33 0.978 83.32 2.6e-04 12.37
4643 COL4A2 13_59 0.976 62.33 1.9e-04 7.76
7786 CATSPER2 15_16 0.973 153.58 4.7e-04 12.68
7960 SERPINF2 17_2 0.964 54.74 1.7e-04 -9.92
9794 GAS2L1 22_10 0.957 26.05 7.9e-05 5.09
12615 EXOC3L2 19_31 0.951 64.16 1.9e-04 -7.88
11318 MIR34AHG 1_6 0.942 23.10 6.9e-05 -4.36
5298 SRP14 15_13 0.940 24.01 7.2e-05 4.54
4401 EIF5A 17_6 0.921 43.88 1.3e-04 -6.32
2376 GALK1 17_42 0.920 35.15 1.0e-04 6.09
11340 PGA3 11_34 0.898 53.44 1.5e-04 -6.35
8045 GLYCTK 3_36 0.893 40.41 1.1e-04 -6.21
9455 PRKAG1 12_31 0.881 39.60 1.1e-04 -6.32
1138 LAT2 7_48 0.875 27.47 7.6e-05 5.28
8815 ANGEL2 1_108 0.855 20.51 5.6e-05 3.91
5637 TPM3 1_75 0.854 31.45 8.5e-05 -6.50
5523 RERE 1_6 0.847 25.39 6.8e-05 4.35
7833 AKTIP 16_29 0.838 22.92 6.1e-05 -4.24
1539 SYNGR1 22_16 0.836 38.59 1.0e-04 5.89
682 EVI5 1_56 0.829 59.90 1.6e-04 -8.00
9179 FAM220A 7_9 0.824 27.46 7.2e-05 -4.69
8537 ORMDL3 17_23 0.817 31.52 8.2e-05 4.73
8337 CD14 5_83 0.808 25.32 6.5e-05 -4.54
#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")
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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
1980 FCGRT 19_34 1.000 238599.30 7.6e-01 17.84
5520 RCN3 19_34 0.000 75878.09 0.0e+00 13.56
8165 CPT1C 19_34 0.000 16519.12 0.0e+00 -6.40
571 SLC6A16 19_34 0.000 4198.39 0.0e+00 -2.89
10492 CTC-301O7.4 19_34 0.000 3959.20 0.0e+00 -1.19
11220 ADM5 19_34 0.000 2502.69 0.0e+00 1.75
6980 ALDH16A1 19_34 0.000 2471.83 0.0e+00 5.37
846 TEAD2 19_34 0.000 2361.64 0.0e+00 0.52
4687 TMEM60 7_49 0.000 1667.57 0.0e+00 1.58
4691 SRPK2 7_65 0.000 1500.51 0.0e+00 2.37
1890 MTFMT 15_30 0.000 1126.36 1.4e-14 3.81
2752 FAM120B 6_112 0.000 1017.80 0.0e+00 2.19
73 KMT2E 7_65 0.000 911.96 0.0e+00 0.83
11489 RP11-325F22.2 7_65 0.000 897.62 0.0e+00 -0.30
1261 SPG21 15_30 0.000 762.94 0.0e+00 -1.47
5519 NOSIP 19_34 0.000 678.44 0.0e+00 -6.15
1989 CD37 19_34 0.000 668.16 0.0e+00 3.39
1541 ASCC2 22_10 0.997 662.80 2.1e-03 3.15
10789 PBX2 6_26 0.000 599.67 0.0e+00 22.53
10825 APOM 6_26 0.000 596.11 0.0e+00 19.97
#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
1980 FCGRT 19_34 1.000 238599.30 0.76000 17.84
1541 ASCC2 22_10 0.997 662.80 0.00210 3.15
6089 FADS1 11_34 0.997 187.54 0.00059 13.66
7786 CATSPER2 15_16 0.973 153.58 0.00047 12.68
19 MAD1L1 7_4 0.319 359.03 0.00036 0.19
7460 ANKRD55 5_33 0.978 83.32 0.00026 12.37
5631 ATP8B2 1_75 0.986 74.88 0.00023 9.37
5712 AFF3 2_58 0.613 106.21 0.00021 10.56
4643 COL4A2 13_59 0.976 62.33 0.00019 7.76
12615 EXOC3L2 19_31 0.951 64.16 0.00019 -7.88
6824 TNFRSF13C 22_17 0.992 60.24 0.00019 7.57
2371 RPS6KB1 17_35 0.774 73.63 0.00018 10.67
7960 SERPINF2 17_2 0.964 54.74 0.00017 -9.92
9073 HIC1 17_3 0.363 139.83 0.00016 12.03
682 EVI5 1_56 0.829 59.90 0.00016 -8.00
11340 PGA3 11_34 0.898 53.44 0.00015 -6.35
4401 EIF5A 17_6 0.921 43.88 0.00013 -6.32
12409 HIST1H3G 6_20 0.677 54.50 0.00012 -1.09
7378 TEX264 3_35 0.564 59.77 0.00011 5.08
9748 TMEM121 14_55 0.419 81.63 0.00011 -8.34
#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
10789 PBX2 6_26 0.000 599.67 0.0e+00 22.53
10825 APOM 6_26 0.000 596.11 0.0e+00 19.97
7712 C2 6_26 0.000 580.67 0.0e+00 -19.92
11047 CLIC1 6_26 0.000 579.66 0.0e+00 19.82
11218 C4B 6_26 0.000 550.82 0.0e+00 -19.79
10808 NEU1 6_26 0.000 575.18 0.0e+00 19.75
11652 C4A 6_26 0.000 538.87 0.0e+00 19.32
1980 FCGRT 19_34 1.000 238599.30 7.6e-01 17.84
10848 TRIM10 6_26 0.000 274.28 2.5e-14 16.58
11374 CYP21A2 6_26 0.000 305.91 0.0e+00 -16.48
4283 TRAF3 14_54 0.013 180.12 7.4e-06 16.30
10781 HLA-DMA 6_27 0.000 220.59 3.3e-09 -15.06
10844 HLA-E 6_26 0.000 530.88 7.3e-10 -14.84
10861 OR2H2 6_23 0.005 111.71 1.8e-06 -14.54
805 FCGR2B 1_79 0.000 500.97 3.6e-09 -14.47
11120 LINC00243 6_26 0.000 143.03 0.0e+00 -13.78
10790 AGER 6_26 0.000 251.94 0.0e+00 -13.69
6089 FADS1 11_34 0.997 187.54 5.9e-04 13.66
5520 RCN3 19_34 0.000 75878.09 0.0e+00 13.56
8659 HSPA6 1_79 0.000 293.67 0.0e+00 -13.34
#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)
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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)
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#proportion of significant z scores
mean(abs(ctwas_gene_res$z) > sig_thresh)
[1] 0.02956287
#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
10789 PBX2 6_26 0.000 599.67 0.0e+00 22.53
10825 APOM 6_26 0.000 596.11 0.0e+00 19.97
7712 C2 6_26 0.000 580.67 0.0e+00 -19.92
11047 CLIC1 6_26 0.000 579.66 0.0e+00 19.82
11218 C4B 6_26 0.000 550.82 0.0e+00 -19.79
10808 NEU1 6_26 0.000 575.18 0.0e+00 19.75
11652 C4A 6_26 0.000 538.87 0.0e+00 19.32
1980 FCGRT 19_34 1.000 238599.30 7.6e-01 17.84
10848 TRIM10 6_26 0.000 274.28 2.5e-14 16.58
11374 CYP21A2 6_26 0.000 305.91 0.0e+00 -16.48
4283 TRAF3 14_54 0.013 180.12 7.4e-06 16.30
10781 HLA-DMA 6_27 0.000 220.59 3.3e-09 -15.06
10844 HLA-E 6_26 0.000 530.88 7.3e-10 -14.84
10861 OR2H2 6_23 0.005 111.71 1.8e-06 -14.54
805 FCGR2B 1_79 0.000 500.97 3.6e-09 -14.47
11120 LINC00243 6_26 0.000 143.03 0.0e+00 -13.78
10790 AGER 6_26 0.000 251.94 0.0e+00 -13.69
6089 FADS1 11_34 0.997 187.54 5.9e-04 13.66
5520 RCN3 19_34 0.000 75878.09 0.0e+00 13.56
8659 HSPA6 1_79 0.000 293.67 0.0e+00 -13.34
ctwas_gene_res_sortz <- ctwas_gene_res[order(-abs(ctwas_gene_res$z)),]
n_plots <- 5
for (region_tag_plot in head(unique(ctwas_gene_res_sortz$region_tag), n_plots)){
ctwas_res_region <- ctwas_res[ctwas_res$region_tag==region_tag_plot,]
start <- min(ctwas_res_region$pos)
end <- max(ctwas_res_region$pos)
ctwas_res_region <- ctwas_res_region[order(ctwas_res_region$pos),]
ctwas_res_region_gene <- ctwas_res_region[ctwas_res_region$type=="gene",]
ctwas_res_region_snp <- ctwas_res_region[ctwas_res_region$type=="SNP",]
#region name
print(paste0("Region: ", region_tag_plot))
#table of genes in region
print(ctwas_res_region_gene[,report_cols])
par(mfrow=c(4,1))
#gene z scores
plot(ctwas_res_region_gene$pos, abs(ctwas_res_region_gene$z), xlab="Position", ylab="abs(gene_z)", xlim=c(start,end),
ylim=c(0,max(sig_thresh, abs(ctwas_res_region_gene$z))),
main=paste0("Region: ", region_tag_plot))
abline(h=sig_thresh,col="red",lty=2)
#significance threshold for SNPs
alpha_snp <- 5*10^(-8)
sig_thresh_snp <- qnorm(1-alpha_snp/2, lower=T)
#snp z scores
plot(ctwas_res_region_snp$pos, abs(ctwas_res_region_snp$z), xlab="Position", ylab="abs(snp_z)",xlim=c(start,end),
ylim=c(0,max(sig_thresh_snp, max(abs(ctwas_res_region_snp$z)))))
abline(h=sig_thresh_snp,col="purple",lty=2)
#gene pips
plot(ctwas_res_region_gene$pos, ctwas_res_region_gene$susie_pip, xlab="Position", ylab="Gene PIP", xlim=c(start,end), ylim=c(0,1))
abline(h=0.8,col="blue",lty=2)
#snp pips
plot(ctwas_res_region_snp$pos, ctwas_res_region_snp$susie_pip, xlab="Position", ylab="SNP PIP", xlim=c(start,end), ylim=c(0,1))
abline(h=0.8,col="blue",lty=2)
}
[1] "Region: 6_26"
genename region_tag susie_pip mu2 PVE z
10855 HLA-G 6_26 0 187.19 0.0e+00 12.67
12599 HCP5B 6_26 0 169.38 0.0e+00 -12.85
10968 HLA-A 6_26 0 234.56 0.0e+00 11.46
10853 HCG9 6_26 0 199.64 0.0e+00 9.31
10851 PPP1R11 6_26 0 51.37 0.0e+00 -1.12
661 ZNRD1 6_26 0 6.38 0.0e+00 -0.14
10850 RNF39 6_26 0 12.47 0.0e+00 2.30
10848 TRIM10 6_26 0 274.28 2.5e-14 16.58
10847 TRIM15 6_26 0 82.44 0.0e+00 6.19
11418 TRIM26 6_26 0 101.32 0.0e+00 -7.15
10845 TRIM39 6_26 0 38.32 0.0e+00 -0.96
11563 RPP21 6_26 0 10.07 0.0e+00 0.77
10844 HLA-E 6_26 0 530.88 7.3e-10 -14.84
10841 MRPS18B 6_26 0 75.93 0.0e+00 2.01
10840 C6orf136 6_26 0 114.21 0.0e+00 -3.68
10839 DHX16 6_26 0 7.54 0.0e+00 0.80
5868 PPP1R18 6_26 0 186.59 0.0e+00 -9.25
4976 NRM 6_26 0 22.75 0.0e+00 2.53
4970 FLOT1 6_26 0 32.83 0.0e+00 -2.66
10230 TUBB 6_26 0 15.71 0.0e+00 0.21
4971 IER3 6_26 0 62.32 0.0e+00 0.06
11120 LINC00243 6_26 0 143.03 0.0e+00 -13.78
10843 DDR1 6_26 0 26.95 0.0e+00 0.98
11052 GTF2H4 6_26 0 5.46 0.0e+00 0.57
4978 VARS2 6_26 0 20.44 0.0e+00 4.24
10838 CCHCR1 6_26 0 85.19 0.0e+00 -6.45
4969 TCF19 6_26 0 73.36 0.0e+00 4.76
10966 HCG27 6_26 0 88.48 0.0e+00 -1.64
10837 POU5F1 6_26 0 145.33 0.0e+00 5.29
10836 HLA-C 6_26 0 69.65 0.0e+00 -5.94
10788 NOTCH4 6_26 0 97.05 0.0e+00 8.77
11439 HLA-B 6_26 0 40.62 0.0e+00 -0.37
12270 XXbac-BPG181B23.7 6_26 0 25.80 0.0e+00 -0.53
10834 MICA 6_26 0 94.57 0.0e+00 -2.59
10833 MICB 6_26 0 20.90 0.0e+00 0.23
10830 LST1 6_26 0 12.29 0.0e+00 2.12
10619 DDX39B 6_26 0 15.38 0.0e+00 4.41
11050 ATP6V1G2 6_26 0 156.09 0.0e+00 -6.87
10831 NFKBIL1 6_26 0 16.67 0.0e+00 -2.93
11282 LTA 6_26 0 36.61 0.0e+00 5.69
11296 LTB 6_26 0 36.30 0.0e+00 5.67
11395 TNF 6_26 0 9.44 0.0e+00 0.93
10829 NCR3 6_26 0 21.58 0.0e+00 -3.84
10828 AIF1 6_26 0 87.63 0.0e+00 -3.04
10827 PRRC2A 6_26 0 76.68 0.0e+00 6.88
10826 BAG6 6_26 0 205.05 0.0e+00 -8.83
10825 APOM 6_26 0 596.11 0.0e+00 19.97
10824 C6orf47 6_26 0 66.79 2.4e-20 -3.01
10822 CSNK2B 6_26 0 82.14 0.0e+00 -8.13
10823 GPANK1 6_26 0 90.67 0.0e+00 9.30
11539 LY6G5B 6_26 0 220.84 0.0e+00 -7.54
10821 LY6G5C 6_26 0 171.76 0.0e+00 -5.70
11639 LY6G6D 6_26 0 187.17 0.0e+00 -6.63
10818 MPIG6B 6_26 0 30.83 0.0e+00 -4.05
10819 LY6G6C 6_26 0 128.47 0.0e+00 -5.91
11048 DDAH2 6_26 0 11.64 0.0e+00 3.97
10817 MSH5 6_26 0 199.80 0.0e+00 -10.02
11047 CLIC1 6_26 0 579.66 0.0e+00 19.82
11327 SAPCD1 6_26 0 52.21 0.0e+00 -0.17
10814 VWA7 6_26 0 33.67 0.0e+00 3.61
10809 C6orf48 6_26 0 71.71 0.0e+00 -3.09
10813 VARS 6_26 0 45.23 0.0e+00 -5.78
10812 LSM2 6_26 0 26.11 0.0e+00 4.29
10811 HSPA1L 6_26 0 73.36 0.0e+00 5.18
10808 NEU1 6_26 0 575.18 0.0e+00 19.75
10807 SLC44A4 6_26 0 182.83 0.0e+00 -10.44
7712 C2 6_26 0 580.67 0.0e+00 -19.92
10805 EHMT2 6_26 0 68.26 0.0e+00 8.37
10802 NELFE 6_26 0 26.24 0.0e+00 4.07
10801 SKIV2L 6_26 0 16.80 0.0e+00 -2.27
10797 STK19 6_26 0 38.19 0.0e+00 5.04
10800 DXO 6_26 0 11.69 0.0e+00 -3.02
11652 C4A 6_26 0 538.87 0.0e+00 19.32
11218 C4B 6_26 0 550.82 0.0e+00 -19.79
11374 CYP21A2 6_26 0 305.91 0.0e+00 -16.48
11193 PPT2 6_26 0 27.32 0.0e+00 -4.81
11043 ATF6B 6_26 0 13.43 0.0e+00 2.66
10795 FKBPL 6_26 0 128.14 0.0e+00 -9.76
10794 PRRT1 6_26 0 92.10 0.0e+00 -9.13
10791 RNF5 6_26 0 65.48 0.0e+00 6.60
11565 EGFL8 6_26 0 28.01 0.0e+00 3.30
10792 AGPAT1 6_26 0 137.97 0.0e+00 10.33
10790 AGER 6_26 0 251.94 0.0e+00 -13.69
10789 PBX2 6_26 0 599.67 0.0e+00 22.53
10608 HLA-DRB5 6_26 0 71.62 0.0e+00 9.50
10325 HLA-DQA1 6_26 0 42.63 0.0e+00 -2.33
11490 HLA-DQA2 6_26 0 74.16 0.0e+00 -4.87
11389 HLA-DQB2 6_26 0 189.47 0.0e+00 -10.57
9260 HLA-DQB1 6_26 0 187.04 0.0e+00 11.07
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
[1] "Region: 19_34"
genename region_tag susie_pip mu2 PVE z
2098 BCAT2 19_34 0 37.38 0.00 1.27
1143 HSD17B14 19_34 0 5.24 0.00 0.02
2100 PLEKHA4 19_34 0 120.15 0.00 2.97
1142 PPP1R15A 19_34 0 115.71 0.00 3.07
1967 NUCB1 19_34 0 17.55 0.00 1.46
1968 DHDH 19_34 0 41.68 0.00 1.33
1146 FTL 19_34 0 211.79 0.00 -2.39
1969 GYS1 19_34 0 40.14 0.00 0.75
9576 RUVBL2 19_34 0 31.91 0.00 -0.87
12166 CTB-60B18.10 19_34 0 65.41 0.00 -0.16
1978 LIN7B 19_34 0 39.51 0.00 -2.79
1976 SNRNP70 19_34 0 449.06 0.00 0.04
11184 C19orf73 19_34 0 176.66 0.00 -0.64
9074 PPFIA3 19_34 0 454.07 0.00 0.43
4200 TRPM4 19_34 0 137.22 0.00 5.18
571 SLC6A16 19_34 0 4198.39 0.00 -2.89
10492 CTC-301O7.4 19_34 0 3959.20 0.00 -1.19
1989 CD37 19_34 0 668.16 0.00 3.39
846 TEAD2 19_34 0 2361.64 0.00 0.52
6980 ALDH16A1 19_34 0 2471.83 0.00 5.37
1980 FCGRT 19_34 1 238599.30 0.76 17.84
5520 RCN3 19_34 0 75878.09 0.00 13.56
5519 NOSIP 19_34 0 678.44 0.00 -6.15
11220 ADM5 19_34 0 2502.69 0.00 1.75
8165 CPT1C 19_34 0 16519.12 0.00 -6.40
3903 TSKS 19_34 0 214.99 0.00 4.80
10357 AP2A1 19_34 0 82.46 0.00 3.01
181 FUZ 19_34 0 23.96 0.00 -1.96
2006 MED25 19_34 0 25.75 0.00 2.49
2001 PTOV1 19_34 0 6.78 0.00 -0.78
381 PNKP 19_34 0 85.04 0.00 -1.50
1998 TBC1D17 19_34 0 31.57 0.00 1.64
8162 ATF5 19_34 0 218.07 0.00 -3.32
10990 NUP62 19_34 0 17.50 0.00 -1.22
6983 SIGLEC11 19_34 0 28.58 0.00 0.80
2015 VRK3 19_34 0 319.43 0.00 1.43
5516 ZNF473 19_34 0 28.53 0.00 1.99
2061 MYH14 19_34 0 247.80 0.00 0.79
4295 NR1H2 19_34 0 21.64 0.00 1.48
4294 NAPSA 19_34 0 46.46 0.00 3.02
6988 EMC10 19_34 0 21.75 0.00 0.73
6989 JOSD2 19_34 0 24.70 0.00 -0.25
10859 ASPDH 19_34 0 50.35 0.00 0.05
2083 CLEC11A 19_34 0 5.61 0.00 0.10
7965 C19orf48 19_34 0 22.99 0.00 0.83
9327 LINC01869 19_34 0 11.52 0.00 -1.65
7966 KLK1 19_34 0 17.70 0.00 -1.54
8152 KLK7 19_34 0 8.90 0.00 0.74
7968 KLK11 19_34 0 20.28 0.00 0.56
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
[1] "Region: 14_54"
genename region_tag susie_pip mu2 PVE z
1242 RCOR1 14_54 0.017 27.14 1.5e-06 -4.74
4283 TRAF3 14_54 0.013 180.12 7.4e-06 16.30
10652 CDC42BPB 14_54 0.009 13.61 3.7e-07 0.05
9776 TNFAIP2 14_54 0.007 6.81 1.5e-07 1.80
1616 EIF5 14_54 0.007 5.37 1.2e-07 0.28
882 MARK3 14_54 0.017 14.22 7.6e-07 -2.20
7687 CKB 14_54 0.013 11.50 4.8e-07 -2.18
7688 TRMT61A 14_54 0.021 15.71 1.0e-06 -2.44
3880 KLC1 14_54 0.031 20.19 2.0e-06 3.22
7689 BAG5 14_54 0.010 9.07 3.0e-07 1.47
11851 APOPT1 14_54 0.010 9.07 3.0e-07 1.47
1191 PPP1R13B 14_54 0.014 12.81 5.7e-07 1.37
6546 TDRD9 14_54 0.009 8.29 2.3e-07 -1.87
6545 C14orf2 14_54 0.030 18.13 1.7e-06 -1.81
7691 ASPG 14_54 0.007 5.47 1.3e-07 0.24
668 KIF26A 14_54 0.018 13.04 7.5e-07 -0.62
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
[1] "Region: 6_27"
genename region_tag susie_pip mu2 PVE z
11559 HLA-DOB 6_27 0.000 9.85 7.0e-12 2.84
10785 TAP2 6_27 0.000 13.60 1.1e-11 -1.93
10782 PSMB8-AS1 6_27 0.000 10.83 2.7e-11 -1.78
10784 PSMB8 6_27 0.000 8.33 4.1e-12 0.39
11540 PSMB9 6_27 0.000 31.78 1.3e-11 -3.99
11596 HLA-DMB 6_27 0.000 86.44 1.9e-08 6.12
10781 HLA-DMA 6_27 0.000 220.59 3.3e-09 -15.06
10780 BRD2 6_27 0.027 94.88 8.1e-06 5.20
10779 HLA-DOA 6_27 0.000 12.24 5.2e-12 -3.30
11209 HLA-DPB1 6_27 0.000 14.97 6.5e-12 -5.30
11360 HLA-DPA1 6_27 0.000 10.98 1.0e-11 -1.45
10778 COL11A2 6_27 0.000 39.29 5.5e-11 -5.57
10775 RXRB 6_27 0.000 23.80 7.0e-11 2.13
10774 HSD17B8 6_27 0.000 8.08 3.8e-12 -1.06
10773 RING1 6_27 0.000 25.83 9.5e-11 2.25
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
[1] "Region: 6_23"
genename region_tag susie_pip mu2 PVE z
10861 OR2H2 6_23 0.005 111.71 1.8e-06 -14.54
10863 GABBR1 6_23 0.002 94.08 5.5e-07 -12.97
10858 ZFP57 6_23 0.002 68.90 3.3e-07 -11.67
10857 HLA-F 6_23 0.001 13.36 6.2e-08 -1.88
10860 MOG 6_23 0.003 59.07 6.2e-07 -8.27
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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
30307 rs61788682 1_69 1.000 39.04 1.2e-04 -5.49
36177 rs61804161 1_79 1.000 448.44 1.4e-03 13.50
36181 rs12145843 1_79 1.000 222.29 7.1e-04 25.22
36188 rs61804205 1_79 1.000 898.37 2.9e-03 -30.68
72900 rs780093 2_16 1.000 222.68 7.1e-04 -16.11
74882 rs17013001 2_21 1.000 34.87 1.1e-04 -5.94
97586 rs62161401 2_66 1.000 70.95 2.3e-04 8.25
152954 rs80351783 3_35 1.000 225.84 7.2e-04 0.85
183113 rs9817452 3_97 1.000 67.16 2.1e-04 8.33
191730 rs9863411 3_114 1.000 65.43 2.1e-04 -8.53
196957 rs3748034 4_4 1.000 79.98 2.5e-04 9.39
277736 rs153429 5_37 1.000 76.85 2.4e-04 -4.54
277751 rs745863029 5_37 1.000 63.97 2.0e-04 -2.08
314270 rs58778501 6_1 1.000 85.83 2.7e-04 -5.80
323592 rs115740542 6_20 1.000 96.25 3.1e-04 -7.06
326714 rs2523581 6_26 1.000 165.53 5.3e-04 -7.25
328434 rs9276685 6_27 1.000 136.60 4.3e-04 -11.75
331732 rs7744080 6_32 1.000 49.95 1.6e-04 7.34
369106 rs60425481 6_104 1.000 764.79 2.4e-03 -5.51
372985 rs139588569 6_112 1.000 9027.76 2.9e-02 -4.62
372987 rs59421548 6_112 1.000 9092.99 2.9e-02 -4.32
373962 rs73041368 7_4 1.000 384.00 1.2e-03 -5.20
391435 rs6583438 7_36 1.000 67.59 2.1e-04 7.15
398517 rs761767938 7_49 1.000 3565.92 1.1e-02 -3.74
398525 rs1544459 7_49 1.000 3501.85 1.1e-02 -3.11
406536 rs763798411 7_65 1.000 11017.62 3.5e-02 3.62
413672 rs125124 7_80 1.000 49.27 1.6e-04 -7.24
467727 rs56114972 8_92 1.000 48.95 1.6e-04 5.92
512810 rs17657502 10_14 1.000 62.54 2.0e-04 5.96
569008 rs12283874 11_36 1.000 59.33 1.9e-04 3.81
580368 rs117304134 11_59 1.000 47.10 1.5e-04 -6.51
583853 rs1176746 11_67 1.000 1337.72 4.2e-03 2.77
583855 rs2307599 11_67 1.000 1336.86 4.2e-03 2.96
635961 rs79490353 13_7 1.000 118.10 3.8e-04 9.78
687606 rs12893029 14_49 1.000 91.89 2.9e-04 -2.30
690666 rs12588969 14_54 1.000 112.51 3.6e-04 -13.24
704099 rs537559727 15_30 1.000 2218.70 7.0e-03 3.09
704108 rs762746560 15_30 1.000 2141.15 6.8e-03 3.19
743610 rs11654694 17_15 1.000 56.55 1.8e-04 -7.82
749524 rs1808192 17_27 1.000 72.23 2.3e-04 9.21
754903 rs113408695 17_39 1.000 48.67 1.5e-04 -7.23
777677 rs150377214 18_35 1.000 74.56 2.4e-04 -8.33
851242 rs11249215 1_17 1.000 57347.86 1.8e-01 -11.47
851248 rs753570588 1_17 1.000 59388.89 1.9e-01 -12.29
1035466 rs62045817 16_54 1.000 53.83 1.7e-04 -7.02
1042083 rs11078597 17_2 1.000 91.76 2.9e-04 11.98
1043277 rs4530175 17_2 1.000 57.05 1.8e-04 7.16
1049391 rs111620634 17_7 1.000 98.81 3.1e-04 7.07
1049464 rs4968200 17_7 1.000 275.81 8.8e-04 10.23
1049519 rs3803800 17_7 1.000 219.50 7.0e-04 -13.97
1074922 rs41523449 19_24 1.000 645.83 2.1e-03 12.98
1074926 rs749726391 19_24 1.000 973.69 3.1e-03 -3.83
1074927 rs12461158 19_24 1.000 938.44 3.0e-03 -4.31
1085523 rs61371437 19_34 1.000 223172.57 7.1e-01 -17.74
1085535 rs374141296 19_34 1.000 225682.83 7.2e-01 -16.57
1108228 rs780018294 22_10 1.000 686.17 2.2e-03 2.17
272258 rs2859493 5_26 0.999 36.43 1.2e-04 6.24
314263 rs4959611 6_1 0.999 60.53 1.9e-04 5.78
427616 rs7012814 8_12 0.999 42.01 1.3e-04 9.08
624567 rs141105880 12_67 0.999 81.67 2.6e-04 -9.92
628487 rs12425627 12_76 0.999 36.63 1.2e-04 -6.14
704103 rs58418704 15_30 0.999 1878.59 6.0e-03 3.38
784467 rs55748813 19_2 0.999 45.46 1.4e-04 -7.13
92554 rs10208803 2_54 0.998 72.97 2.3e-04 7.79
97476 rs12622400 2_66 0.998 42.28 1.3e-04 5.70
252647 rs7659414 4_114 0.998 47.10 1.5e-04 -7.21
292126 rs35552666 5_66 0.998 32.07 1.0e-04 -5.80
485007 rs4745108 9_33 0.998 30.96 9.8e-05 -5.45
586551 rs1945396 11_75 0.998 34.81 1.1e-04 5.84
643023 rs7997446 13_21 0.998 32.98 1.0e-04 6.01
704106 rs11858985 15_30 0.997 2104.44 6.7e-03 2.96
706795 rs56357772 15_36 0.997 60.49 1.9e-04 -9.14
221843 rs12507099 4_53 0.996 29.80 9.4e-05 -5.41
234118 rs138204164 4_77 0.996 60.62 1.9e-04 -7.95
324730 rs1233385 6_23 0.996 126.42 4.0e-04 -14.94
535122 rs10887917 10_57 0.996 48.53 1.5e-04 7.04
334869 rs6458803 6_39 0.995 33.14 1.0e-04 5.71
628408 rs2229840 12_75 0.995 30.52 9.6e-05 -5.48
705471 rs35485240 15_33 0.995 68.63 2.2e-04 -8.57
49923 rs3813977 1_105 0.994 33.82 1.1e-04 5.48
84697 rs13012253 2_39 0.994 29.13 9.2e-05 -5.37
132402 rs1834748 2_135 0.994 36.36 1.1e-04 6.44
748695 rs8072356 17_26 0.993 31.08 9.8e-05 5.11
833769 rs12166267 22_6 0.993 32.44 1.0e-04 5.49
192609 rs79692229 3_116 0.992 39.92 1.3e-04 6.67
1049849 rs78378222 17_7 0.991 64.94 2.0e-04 -4.72
229187 rs144812644 4_68 0.990 29.32 9.2e-05 -6.53
687623 rs2239651 14_49 0.990 50.89 1.6e-04 -7.09
361346 rs6924387 6_90 0.989 88.29 2.8e-04 9.78
514376 rs148678804 10_16 0.989 27.84 8.7e-05 4.86
258828 rs56023411 5_2 0.988 37.45 1.2e-04 6.31
323571 rs72834643 6_20 0.988 33.01 1.0e-04 -3.68
691596 rs61310292 14_56 0.988 48.77 1.5e-04 -7.41
724923 rs8061729 16_24 0.988 39.63 1.2e-04 5.12
759555 rs9954032 18_1 0.988 31.07 9.7e-05 -5.43
585514 rs666741 11_71 0.984 63.89 2.0e-04 -8.52
668833 rs8011368 14_10 0.984 27.29 8.5e-05 5.02
622149 rs2583223 12_62 0.983 34.31 1.1e-04 -5.70
776411 rs12960077 18_32 0.983 34.38 1.1e-04 -5.85
791006 rs71332143 19_15 0.983 28.80 9.0e-05 -5.44
328296 rs112357706 6_27 0.982 38.71 1.2e-04 5.80
426622 rs11775663 8_10 0.982 27.90 8.7e-05 -5.28
463108 rs2720659 8_84 0.981 33.61 1.0e-04 -5.88
765151 rs35796589 18_10 0.981 26.13 8.1e-05 -4.64
544709 rs11199973 10_75 0.980 31.56 9.8e-05 -5.52
124127 rs231811 2_120 0.979 53.18 1.7e-04 7.75
324893 rs2246856 6_23 0.978 60.97 1.9e-04 -5.56
611947 rs2137537 12_44 0.978 29.96 9.3e-05 -5.45
498941 rs2812398 9_58 0.975 31.31 9.7e-05 5.52
399921 rs3839804 7_51 0.972 29.41 9.1e-05 -5.44
568949 rs11227230 11_36 0.971 53.48 1.6e-04 -5.29
673780 rs2883893 14_20 0.969 28.06 8.6e-05 -5.85
303412 rs12189018 5_87 0.968 25.91 8.0e-05 4.88
777717 rs4940573 18_35 0.967 119.33 3.7e-04 10.74
141507 rs17776482 3_9 0.966 27.03 8.3e-05 -5.33
280232 rs6452453 5_43 0.962 27.25 8.3e-05 -5.25
687714 rs2069987 14_49 0.962 31.69 9.7e-05 -5.58
196958 rs3752442 4_4 0.958 48.43 1.5e-04 -9.89
428288 rs7821812 8_14 0.958 95.13 2.9e-04 -11.91
126117 rs62203749 2_124 0.957 25.96 7.9e-05 -4.49
426828 rs12543422 8_10 0.955 25.08 7.6e-05 4.59
409769 rs38913 7_71 0.953 27.76 8.4e-05 5.22
694402 rs7497631 15_7 0.951 24.65 7.4e-05 -4.59
74278 rs115472871 2_20 0.947 24.92 7.5e-05 -4.83
607906 rs7397189 12_36 0.947 26.60 8.0e-05 -4.79
702721 rs340029 15_27 0.942 58.77 1.8e-04 7.80
314286 rs6942338 6_1 0.938 84.13 2.5e-04 10.28
588625 rs7932045 11_80 0.937 30.53 9.1e-05 7.01
324077 rs187257713 6_21 0.935 24.95 7.4e-05 -3.89
785574 rs67868323 19_4 0.935 61.57 1.8e-04 8.04
482594 rs11557154 9_27 0.932 35.73 1.1e-04 -5.67
30305 rs56894897 1_69 0.931 26.22 7.8e-05 -3.72
504066 rs495828 9_70 0.930 36.85 1.1e-04 5.32
831294 rs12626883 21_24 0.929 24.62 7.3e-05 -4.68
635963 rs7989654 13_7 0.928 63.62 1.9e-04 5.76
141725 rs56395424 3_9 0.927 32.09 9.5e-05 -4.58
512825 rs2497836 10_14 0.926 40.34 1.2e-04 -3.55
382657 rs10228771 7_21 0.925 24.04 7.1e-05 -4.46
196962 rs1203107 4_4 0.924 70.79 2.1e-04 8.52
370233 rs9365555 6_106 0.923 24.24 7.1e-05 -4.61
370336 rs766167 6_106 0.923 25.14 7.4e-05 -4.85
570515 rs11235597 11_41 0.919 24.52 7.2e-05 -4.49
687418 rs12588988 14_47 0.916 23.92 7.0e-05 4.62
664887 rs35477689 14_3 0.914 40.37 1.2e-04 -6.86
680931 rs61987084 14_34 0.914 28.64 8.3e-05 -5.11
470904 rs10120959 9_4 0.913 23.77 6.9e-05 -4.51
318363 rs45449792 6_10 0.912 23.61 6.8e-05 4.49
379985 rs111683935 7_17 0.912 31.54 9.1e-05 -5.58
799993 rs670795 19_37 0.911 47.32 1.4e-04 -6.95
74919 rs13388394 2_21 0.910 25.93 7.5e-05 -5.05
33328 rs12124727 1_73 0.909 25.93 7.5e-05 3.48
624492 rs653178 12_67 0.909 55.15 1.6e-04 -8.35
508169 rs1972409 10_7 0.906 34.37 9.9e-05 6.24
152967 rs71329026 3_35 0.902 214.52 6.1e-04 3.37
10219 rs2045791 1_23 0.899 23.64 6.7e-05 -4.39
323551 rs140264349 6_20 0.899 31.70 9.1e-05 -4.78
63710 rs4335411 1_131 0.896 23.88 6.8e-05 -4.41
213741 rs768294452 4_39 0.896 23.64 6.7e-05 3.88
1080145 rs148933445 19_31 0.896 32.42 9.2e-05 -5.54
398521 rs11972122 7_49 0.895 3288.10 9.3e-03 -3.63
814251 rs74178731 20_29 0.894 28.71 8.2e-05 5.34
406542 rs13230660 7_65 0.893 10947.58 3.1e-02 4.37
427126 rs7833103 8_11 0.890 37.43 1.1e-04 7.32
665219 rs12588750 14_3 0.887 35.74 1.0e-04 -5.86
594589 rs10743892 12_10 0.886 40.44 1.1e-04 -6.22
504444 rs7043538 9_71 0.884 25.10 7.0e-05 -4.64
875534 rs13063578 3_33 0.883 51.26 1.4e-04 -6.79
349875 rs2388334 6_67 0.881 36.18 1.0e-04 -5.94
588597 rs6590334 11_80 0.879 38.02 1.1e-04 7.28
593395 rs10734885 12_7 0.879 26.91 7.5e-05 -4.76
99867 rs2422391 2_69 0.877 28.44 7.9e-05 -5.00
499295 rs2418317 9_59 0.877 30.58 8.5e-05 -5.36
229157 rs17032996 4_68 0.875 32.06 8.9e-05 6.69
8377 rs2491141 1_20 0.874 24.94 6.9e-05 4.71
743633 rs3751985 17_15 0.873 444.65 1.2e-03 25.61
148638 rs116823501 3_24 0.870 23.87 6.6e-05 3.20
536163 rs11187129 10_59 0.870 30.73 8.5e-05 3.80
504398 rs56406717 9_70 0.868 25.55 7.0e-05 -4.87
1048635 rs148272371 17_6 0.865 30.02 8.2e-05 -5.05
593495 rs4883268 12_7 0.855 28.86 7.8e-05 -4.98
224510 rs114646961 4_59 0.853 24.42 6.6e-05 4.50
188991 rs2141598 3_109 0.849 24.09 6.5e-05 4.45
47473 rs112840522 1_99 0.847 23.87 6.4e-05 -4.19
553723 rs360130 11_8 0.846 44.42 1.2e-04 -5.27
107293 rs60882035 2_85 0.844 33.23 8.9e-05 -6.19
133492 rs34013402 2_137 0.844 35.99 9.6e-05 -6.00
496771 rs10733564 9_54 0.844 24.79 6.6e-05 4.45
511869 rs10906857 10_13 0.838 23.51 6.3e-05 -4.31
752152 rs1040261 17_33 0.838 29.76 7.9e-05 -5.32
832169 rs62222326 22_3 0.838 28.10 7.5e-05 -5.05
1043201 rs11651121 17_2 0.835 34.06 9.0e-05 5.22
294770 rs71583081 5_71 0.834 23.54 6.2e-05 -4.32
54094 rs61830291 1_112 0.831 48.37 1.3e-04 7.06
373875 rs79206451 7_3 0.831 26.31 6.9e-05 -4.55
661882 rs17381234 13_57 0.831 25.83 6.8e-05 4.74
178167 rs9862179 3_86 0.829 24.79 6.5e-05 -4.46
1006273 rs5742915 15_35 0.826 30.72 8.1e-05 -5.11
47087 rs74213209 1_98 0.825 34.77 9.1e-05 -5.83
275582 rs28499105 5_31 0.825 31.39 8.2e-05 5.33
659838 rs35832914 13_52 0.825 28.67 7.5e-05 5.14
365861 rs10872678 6_99 0.823 32.87 8.6e-05 5.60
714666 rs10902585 15_49 0.823 24.90 6.5e-05 4.58
205674 rs10007850 4_22 0.820 144.50 3.8e-04 2.08
374217 rs12671734 7_5 0.817 26.10 6.8e-05 4.50
569197 rs574546203 11_37 0.817 25.76 6.7e-05 4.62
815813 rs140571612 20_32 0.817 24.05 6.2e-05 -4.33
488889 rs930340 9_41 0.812 60.78 1.6e-04 -8.21
568944 rs59286748 11_36 0.810 42.03 1.1e-04 -5.94
196940 rs115019205 4_4 0.809 26.56 6.8e-05 4.69
361689 rs112020444 6_92 0.807 42.39 1.1e-04 6.57
96063 rs6711659 2_63 0.806 25.69 6.6e-05 4.59
329671 rs9348980 6_29 0.803 34.75 8.9e-05 5.33
328027 rs9273504 6_26 0.802 661.33 1.7e-03 -24.74
#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")
Version | Author | Date |
---|---|---|
dfd2b5f | wesleycrouse | 2021-09-07 |
#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
1085535 rs374141296 19_34 1 225682.8 0.72 -16.57
1085532 rs113176985 19_34 0 223722.4 0.00 -17.85
1085539 rs2946865 19_34 0 223313.9 0.00 -17.84
1085525 rs35295508 19_34 0 223273.0 0.00 -17.73
1085523 rs61371437 19_34 1 223172.6 0.71 -17.74
1085530 rs73056069 19_34 0 222510.0 0.00 -17.77
1085514 rs756628 19_34 0 222259.8 0.00 -17.67
1085513 rs739349 19_34 0 222258.9 0.00 -17.67
1085527 rs2878354 19_34 0 221933.9 0.00 -17.74
1085510 rs739347 19_34 0 221821.3 0.00 -17.75
1085511 rs2073614 19_34 0 221548.1 0.00 -17.75
1085516 rs2077300 19_34 0 220958.7 0.00 -17.69
1085520 rs73056059 19_34 0 220560.7 0.00 -17.73
1085506 rs4802613 19_34 0 220557.6 0.00 -17.66
1085540 rs60815603 19_34 0 219522.9 0.00 -17.84
1085543 rs1316885 19_34 0 219092.2 0.00 -17.95
1085548 rs2946863 19_34 0 218715.9 0.00 -18.01
1085541 rs35443645 19_34 0 218418.6 0.00 -18.08
1085545 rs60746284 19_34 0 218190.6 0.00 -17.92
1085504 rs10403394 19_34 0 217555.7 0.00 -17.58
1085505 rs17555056 19_34 0 217409.8 0.00 -17.63
1085521 rs73056062 19_34 0 214812.7 0.00 -17.00
1085551 rs553431297 19_34 0 212066.3 0.00 -17.07
1085534 rs112283514 19_34 0 211330.4 0.00 -16.37
1085536 rs11270139 19_34 0 210057.6 0.00 -16.89
1085501 rs10421294 19_34 0 196549.9 0.00 -16.76
1085500 rs8108175 19_34 0 196523.0 0.00 -16.77
1085493 rs59192944 19_34 0 196147.8 0.00 -16.75
1085499 rs1858742 19_34 0 196079.7 0.00 -16.77
1085490 rs55991145 19_34 0 195983.9 0.00 -16.79
1085485 rs3786567 19_34 0 195906.3 0.00 -16.78
1085481 rs2271952 19_34 0 195828.0 0.00 -16.79
1085484 rs4801801 19_34 0 195804.2 0.00 -16.79
1085480 rs2271953 19_34 0 195588.3 0.00 -16.82
1085482 rs2271951 19_34 0 195580.2 0.00 -16.82
1085471 rs60365978 19_34 0 195433.9 0.00 -16.84
1085477 rs4802612 19_34 0 194660.8 0.00 -16.81
1085487 rs2517977 19_34 0 194391.7 0.00 -16.80
1085474 rs55893003 19_34 0 194140.1 0.00 -16.85
1085466 rs55992104 19_34 0 189480.2 0.00 -15.98
1085460 rs60403475 19_34 0 189458.3 0.00 -15.96
1085463 rs4352151 19_34 0 189398.2 0.00 -15.99
1085457 rs11878448 19_34 0 189261.9 0.00 -15.98
1085451 rs9653100 19_34 0 189223.3 0.00 -15.97
1085447 rs4802611 19_34 0 189101.6 0.00 -15.98
1085439 rs7251338 19_34 0 188814.5 0.00 -15.99
1085438 rs59269605 19_34 0 188792.7 0.00 -16.00
1085459 rs1042120 19_34 0 188304.7 0.00 -16.00
1085455 rs113220577 19_34 0 188138.8 0.00 -16.00
1085449 rs9653118 19_34 0 187842.2 0.00 -15.98
#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
1085535 rs374141296 19_34 1.000 225682.83 0.72000 -16.57
1085523 rs61371437 19_34 1.000 223172.57 0.71000 -17.74
851248 rs753570588 1_17 1.000 59388.89 0.19000 -12.29
851242 rs11249215 1_17 1.000 57347.86 0.18000 -11.47
851267 rs11249219 1_17 0.792 57208.95 0.14000 -13.55
851258 rs12407074 1_17 0.306 57231.49 0.05600 -13.56
851256 rs7555518 1_17 0.211 57231.19 0.03800 -13.55
406536 rs763798411 7_65 1.000 11017.62 0.03500 3.62
851264 rs7513156 1_17 0.183 57230.38 0.03300 -13.55
406542 rs13230660 7_65 0.893 10947.58 0.03100 4.37
372985 rs139588569 6_112 1.000 9027.76 0.02900 -4.62
372987 rs59421548 6_112 1.000 9092.99 0.02900 -4.32
406547 rs4997569 7_65 0.752 10972.10 0.02600 4.29
851265 rs10903121 1_17 0.130 57229.71 0.02400 -13.55
851261 rs7550635 1_17 0.124 57230.58 0.02300 -13.55
406554 rs6952534 7_65 0.571 10946.79 0.02000 4.44
851263 rs7542123 1_17 0.108 57230.09 0.02000 -13.55
851260 rs7550552 1_17 0.104 57230.23 0.01900 -13.55
406539 rs10274607 7_65 0.376 10962.76 0.01300 4.32
398517 rs761767938 7_49 1.000 3565.92 0.01100 -3.74
398525 rs1544459 7_49 1.000 3501.85 0.01100 -3.11
398521 rs11972122 7_49 0.895 3288.10 0.00930 -3.63
704099 rs537559727 15_30 1.000 2218.70 0.00700 3.09
704108 rs762746560 15_30 1.000 2141.15 0.00680 3.19
704106 rs11858985 15_30 0.997 2104.44 0.00670 2.96
704103 rs58418704 15_30 0.999 1878.59 0.00600 3.38
406553 rs4730069 7_65 0.147 10939.21 0.00510 4.45
398464 rs9640663 7_49 0.635 2104.85 0.00420 -5.07
583853 rs1176746 11_67 1.000 1337.72 0.00420 2.77
583855 rs2307599 11_67 1.000 1336.86 0.00420 2.96
1074926 rs749726391 19_24 1.000 973.69 0.00310 -3.83
1074927 rs12461158 19_24 1.000 938.44 0.00300 -4.31
36188 rs61804205 1_79 1.000 898.37 0.00290 -30.68
369106 rs60425481 6_104 1.000 764.79 0.00240 -5.51
398460 rs2868787 7_49 0.365 2103.94 0.00240 -5.02
1108228 rs780018294 22_10 1.000 686.17 0.00220 2.17
1074922 rs41523449 19_24 1.000 645.83 0.00210 12.98
328027 rs9273504 6_26 0.802 661.33 0.00170 -24.74
406546 rs10242713 7_65 0.046 10908.93 0.00160 4.52
36177 rs61804161 1_79 1.000 448.44 0.00140 13.50
373962 rs73041368 7_4 1.000 384.00 0.00120 -5.20
743633 rs3751985 17_15 0.873 444.65 0.00120 25.61
398522 rs11406602 7_49 0.105 3282.18 0.00110 -3.52
1074961 rs45512696 19_24 0.391 782.43 0.00097 14.72
1049464 rs4968200 17_7 1.000 275.81 0.00088 10.23
851266 rs11249218 1_17 0.005 57215.21 0.00083 -13.53
1074963 rs58895965 19_24 0.321 781.84 0.00080 14.72
36137 rs115032752 1_79 0.664 347.51 0.00073 -23.12
152954 rs80351783 3_35 1.000 225.84 0.00072 0.85
36181 rs12145843 1_79 1.000 222.29 0.00071 25.22
#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
36188 rs61804205 1_79 1.000 898.37 2.9e-03 -30.68
36164 rs189026820 1_79 0.000 766.03 0.0e+00 -27.18
36184 rs74816838 1_79 0.000 648.93 0.0e+00 -26.47
36143 rs7518087 1_79 0.000 499.24 8.8e-19 -25.90
1085645 rs36013629 19_34 0.000 44412.13 0.0e+00 -25.88
1085620 rs10419198 19_34 0.000 45152.82 0.0e+00 -25.81
743633 rs3751985 17_15 0.873 444.65 1.2e-03 25.61
743631 rs3794776 17_15 0.129 454.46 1.9e-04 25.37
36181 rs12145843 1_79 1.000 222.29 7.1e-04 25.22
743628 rs16961828 17_15 0.002 437.18 2.3e-06 25.04
328027 rs9273504 6_26 0.802 661.33 1.7e-03 -24.74
328068 rs9274465 6_26 0.198 658.31 4.1e-04 -24.68
328023 rs9273455 6_26 0.000 647.70 7.7e-07 -24.49
1085703 rs111476047 19_34 0.000 41765.55 0.0e+00 -23.78
328048 rs17613606 6_26 0.000 638.74 6.9e-12 -23.70
36159 rs61801830 1_79 0.000 218.92 2.1e-09 -23.65
36137 rs115032752 1_79 0.664 347.51 7.3e-04 -23.12
36134 rs146188788 1_79 0.336 345.32 3.7e-04 -23.07
327155 rs3130281 6_26 0.000 613.56 6.5e-19 -22.62
327160 rs3131297 6_26 0.000 613.28 6.5e-19 -22.61
1085553 rs2335534 19_34 0.000 113531.91 0.0e+00 -22.38
1085577 rs10469298 19_34 0.000 68996.57 0.0e+00 -22.06
1085590 rs1132990 19_34 0.000 69339.44 0.0e+00 -21.97
36133 rs9427059 1_79 0.000 305.04 0.0e+00 -21.81
1085771 rs2379087 19_34 0.000 33853.67 0.0e+00 -21.56
1085765 rs111310942 19_34 0.000 33984.92 0.0e+00 -21.51
1085820 rs11878568 19_34 0.000 33375.12 0.0e+00 -21.48
1085775 rs10416310 19_34 0.000 33650.27 0.0e+00 -21.47
1085811 rs7249925 19_34 0.000 33340.51 0.0e+00 -21.46
1085781 rs3760708 19_34 0.000 33193.88 0.0e+00 -21.45
1085763 rs7254718 19_34 0.000 34361.90 0.0e+00 -21.44
1085834 rs3745475 19_34 0.000 33135.09 0.0e+00 -21.44
1085840 rs10417980 19_34 0.000 35515.49 0.0e+00 -21.44
1085827 rs10414643 19_34 0.000 34873.78 0.0e+00 -21.43
1085785 rs10421333 19_34 0.000 32940.07 0.0e+00 -21.41
1085792 rs2890072 19_34 0.000 32850.99 0.0e+00 -21.41
1085798 rs10426059 19_34 0.000 32229.18 0.0e+00 -21.41
1085801 rs112727702 19_34 0.000 32682.14 0.0e+00 -21.41
1085795 rs8113357 19_34 0.000 32724.15 0.0e+00 -21.40
1085803 rs10406941 19_34 0.000 32768.40 0.0e+00 -21.40
1085806 rs2288921 19_34 0.000 32728.54 0.0e+00 -21.40
1085842 rs2116922 19_34 0.000 35658.38 0.0e+00 -21.40
1085791 rs2379088 19_34 0.000 32806.85 0.0e+00 -21.39
1085821 rs11083979 19_34 0.000 32747.15 0.0e+00 -21.39
1085823 rs7251295 19_34 0.000 32722.26 0.0e+00 -21.39
1085825 rs7251877 19_34 0.000 32561.71 0.0e+00 -21.39
1085800 rs56873913 19_34 0.000 32752.66 0.0e+00 -21.38
1085822 rs10404887 19_34 0.000 32699.08 0.0e+00 -21.38
1085828 rs16981329 19_34 0.000 32643.09 0.0e+00 -21.37
1085744 rs10412446 19_34 0.000 34391.79 0.0e+00 -21.36
#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.8]
#number of genes for gene set enrichment
length(genes)
[1] 29
if (length(genes)>0){
GO_enrichment <- enrichr(genes, dbs)
for (db in dbs){
print(db)
df <- 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"
[1] Term Overlap Adjusted.P.value Genes
<0 rows> (or 0-length row.names)
[1] "GO_Cellular_Component_2021"
[1] Term Overlap Adjusted.P.value Genes
<0 rows> (or 0-length row.names)
[1] "GO_Molecular_Function_2021"
[1] Term Overlap Adjusted.P.value Genes
<0 rows> (or 0-length row.names)
CHMP7 gene(s) from the input list not found in DisGeNET CURATEDANGEL2 gene(s) from the input list not found in DisGeNET CURATEDFAM220A gene(s) from the input list not found in DisGeNET CURATEDPRKAG1 gene(s) from the input list not found in DisGeNET CURATEDLAT2 gene(s) from the input list not found in DisGeNET CURATEDPGA3 gene(s) from the input list not found in DisGeNET CURATEDSRP14 gene(s) from the input list not found in DisGeNET CURATEDMIR34AHG gene(s) from the input list not found in DisGeNET CURATEDATP8B2 gene(s) from the input list not found in DisGeNET CURATEDEVI5 gene(s) from the input list not found in DisGeNET CURATEDFCGRT gene(s) from the input list not found in DisGeNET CURATEDASCC2 gene(s) from the input list not found in DisGeNET CURATED
Description FDR Ratio BgRatio
48 Renal Cell Dysplasia 0.01496123 1/17 1/9703
57 D-Glyceric aciduria 0.01496123 1/17 1/9703
65 Anhydramnios 0.01496123 1/17 1/9703
72 D-glycericacidemia 0.01496123 1/17 1/9703
76 Nemaline myopathy 1 0.01496123 1/17 1/9703
82 CAP MYOPATHY, TPM3-RELATED (disorder) 0.01496123 1/17 1/9703
84 ALPHA-2-PLASMIN INHIBITOR DEFICIENCY 0.01496123 1/17 1/9703
89 IMMUNODEFICIENCY, COMMON VARIABLE, 4 0.01496123 1/17 1/9703
91 PORENCEPHALY 2 0.01496123 1/17 1/9703
93 Anti-plasmin deficiency, congenital 0.01496123 1/17 1/9703
******************************************
* *
* Welcome to WebGestaltR ! *
* *
******************************************
Loading the functional categories...
Loading the ID list...
Loading the reference list...
Performing the enrichment analysis...
Warning in oraEnrichment(interestGeneList, referenceGeneList, geneSet,
minNum = minNum, : No significant gene set is identified based on FDR 0.05!
NULL
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] WebGestaltR_0.4.4 disgenet2r_0.99.2 enrichR_3.0 cowplot_1.0.0
[5] ggplot2_3.3.3
loaded via a namespace (and not attached):
[1] bitops_1.0-6 matrixStats_0.57.0
[3] fs_1.3.1 bit64_4.0.5
[5] doParallel_1.0.16 progress_1.2.2
[7] httr_1.4.1 rprojroot_2.0.2
[9] GenomeInfoDb_1.20.0 doRNG_1.8.2
[11] tools_3.6.1 utf8_1.2.1
[13] R6_2.5.0 DBI_1.1.1
[15] BiocGenerics_0.30.0 colorspace_1.4-1
[17] withr_2.4.1 tidyselect_1.1.0
[19] prettyunits_1.0.2 bit_4.0.4
[21] curl_3.3 compiler_3.6.1
[23] git2r_0.26.1 Biobase_2.44.0
[25] DelayedArray_0.10.0 rtracklayer_1.44.0
[27] labeling_0.3 scales_1.1.0
[29] readr_1.4.0 apcluster_1.4.8
[31] stringr_1.4.0 digest_0.6.20
[33] Rsamtools_2.0.0 svglite_1.2.2
[35] rmarkdown_1.13 XVector_0.24.0
[37] pkgconfig_2.0.3 htmltools_0.3.6
[39] fastmap_1.1.0 BSgenome_1.52.0
[41] rlang_0.4.11 RSQLite_2.2.7
[43] generics_0.0.2 farver_2.1.0
[45] jsonlite_1.6 BiocParallel_1.18.0
[47] dplyr_1.0.7 VariantAnnotation_1.30.1
[49] RCurl_1.98-1.1 magrittr_2.0.1
[51] GenomeInfoDbData_1.2.1 Matrix_1.2-18
[53] Rcpp_1.0.6 munsell_0.5.0
[55] S4Vectors_0.22.1 fansi_0.5.0
[57] gdtools_0.1.9 lifecycle_1.0.0
[59] stringi_1.4.3 whisker_0.3-2
[61] yaml_2.2.0 SummarizedExperiment_1.14.1
[63] zlibbioc_1.30.0 plyr_1.8.4
[65] grid_3.6.1 blob_1.2.1
[67] parallel_3.6.1 promises_1.0.1
[69] crayon_1.4.1 lattice_0.20-38
[71] Biostrings_2.52.0 GenomicFeatures_1.36.3
[73] hms_1.1.0 knitr_1.23
[75] pillar_1.6.1 igraph_1.2.4.1
[77] GenomicRanges_1.36.0 rjson_0.2.20
[79] rngtools_1.5 codetools_0.2-16
[81] reshape2_1.4.3 biomaRt_2.40.1
[83] stats4_3.6.1 XML_3.98-1.20
[85] glue_1.4.2 evaluate_0.14
[87] data.table_1.14.0 foreach_1.5.1
[89] vctrs_0.3.8 httpuv_1.5.1
[91] gtable_0.3.0 purrr_0.3.4
[93] assertthat_0.2.1 cachem_1.0.5
[95] xfun_0.8 later_0.8.0
[97] tibble_3.1.2 iterators_1.0.13
[99] GenomicAlignments_1.20.1 AnnotationDbi_1.46.0
[101] memoise_2.0.0 IRanges_2.18.1
[103] workflowr_1.6.2 ellipsis_0.3.2