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

Overview

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])

Weight QC

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

#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)

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)

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

Genes with highest PIPs

#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

Genes 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")

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 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
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 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
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

Comparing 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)

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

Locus plots for genes and SNPs

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 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
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

SNPs 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")

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 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
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 largest z scores

#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

Gene 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.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