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-30870_irnt_Whole_Blood.Rmd
    Modified:   analysis/ukb-d-30880_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30880_irnt_Whole_Blood.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-30880_irnt_Whole_Blood.Rmd) and HTML (docs/ukb-d-30880_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 Urate (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-30880_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.008022976 0.000204890 
#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 
20.81243 18.20088 
#report sample size
print(sample_size)
[1] 343836
#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.005388083 0.094329551 
#compare sum(PIP*mu2/sample_size) with above PVE calculation
c(sum(ctwas_gene_res$PVE),sum(ctwas_snp_res$PVE))
[1] 0.01825411 0.71407143

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
10765  ZDHHC18       1_18     1.000 142.25 4.1e-04  12.20
9813      MUC1       1_77     1.000 246.05 7.2e-04 -18.55
3752    KCNK17       6_30     0.971  39.18 1.1e-04  -6.59
10816     MBD5       2_88     0.969  81.94 2.3e-04   9.34
6598      NRG1       8_31     0.966  71.88 2.0e-04  -8.59
5012    TRIM29      11_72     0.941  47.36 1.3e-04   7.00
8523    ZNF217      20_31     0.932  25.40 6.9e-05  -4.69
5874      ANO7      2_143     0.913  37.89 1.0e-04  -6.01
6726    ZNF276      16_54     0.894  31.37 8.2e-05  -5.44
9820     TLCD2       17_2     0.891  25.23 6.5e-05  -4.79
8177     THBS3       1_77     0.890 216.19 5.6e-04  16.74
4484   CSNK1G2       19_2     0.859  27.97 7.0e-05   5.00
10355 SLC39A10      2_116     0.846  27.35 6.7e-05  -4.94
4505     AGAP3       7_94     0.823  25.21 6.0e-05  -4.56
9325    ZNF816      19_37     0.803  33.38 7.8e-05  -5.73
920     STXBP2       19_7     0.798  23.95 5.6e-05   4.38
1916    TRIM35       8_27     0.772  24.88 5.6e-05   4.74
2540   NECTIN1      11_72     0.753  22.58 4.9e-05  -4.47
5564    ATP1B1       1_82     0.751  24.21 5.3e-05   4.27
12181    EGLN2      19_30     0.746  24.64 5.3e-05  -4.12

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
6501     PDIA4       7_92     0.000 7907.58 0.0e+00   2.15
168      SPRTN      1_118     0.000 5408.85 1.6e-10   3.40
3138     EXOC8      1_118     0.000 3884.27 0.0e+00   2.52
4687    TMEM60       7_49     0.000 2850.57 0.0e+00   2.04
9157   ZNF518B       4_11     0.000 1433.04 0.0e+00  58.12
2475    SLC2A9       4_11     0.000 1391.98 0.0e+00  53.26
10421   ZNF786       7_92     0.000 1280.98 0.0e+00  -1.75
8305    ZNF282       7_92     0.000 1229.89 0.0e+00  -1.63
10367   ZNF398       7_92     0.000  948.52 0.0e+00  -1.27
3140     TSNAX      1_118     0.000  758.59 0.0e+00   1.97
11094     APTR       7_49     0.000  494.47 0.0e+00  -2.79
7145     DISC1      1_118     0.000  459.53 0.0e+00  -1.78
98       PHTF2       7_49     0.000  415.54 0.0e+00   0.42
10882   ZNF425       7_92     0.000  381.88 0.0e+00  -0.60
2737    TRIM38       6_20     0.005  315.61 4.5e-06 -23.44
9813      MUC1       1_77     1.000  246.05 7.2e-04 -18.55
10864    SPDYC      11_36     0.001  230.89 9.2e-07  -9.34
8177     THBS3       1_77     0.890  216.19 5.6e-04  16.74
2953     NRBP1       2_16     0.021  215.65 1.3e-05  15.35
11023    SIPA1      11_36     0.000  214.21 8.6e-11  12.60

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
9813          MUC1       1_77     1.000 246.05 7.2e-04 -18.55
8177         THBS3       1_77     0.890 216.19 5.6e-04  16.74
10765      ZDHHC18       1_18     1.000 142.25 4.1e-04  12.20
10816         MBD5       2_88     0.969  81.94 2.3e-04   9.34
6598          NRG1       8_31     0.966  71.88 2.0e-04  -8.59
5012        TRIM29      11_72     0.941  47.36 1.3e-04   7.00
3752        KCNK17       6_30     0.971  39.18 1.1e-04  -6.59
6089         FADS1      11_34     0.679  52.66 1.0e-04   7.27
5874          ANO7      2_143     0.913  37.89 1.0e-04  -6.01
171         UQCRC1       3_34     0.687  43.25 8.6e-05  -6.52
3544        FAM35A      10_57     0.740  39.04 8.4e-05  -6.62
6726        ZNF276      16_54     0.894  31.37 8.2e-05  -5.44
9325        ZNF816      19_37     0.803  33.38 7.8e-05  -5.73
12389 RP11-686O6.2      2_120     0.359  71.39 7.5e-05   8.04
10021      ZKSCAN4       6_22     0.155 154.44 7.0e-05  12.61
4484       CSNK1G2       19_2     0.859  27.97 7.0e-05   5.00
8523        ZNF217      20_31     0.932  25.40 6.9e-05  -4.69
10355     SLC39A10      2_116     0.846  27.35 6.7e-05  -4.94
4126        QRICH2      17_42     0.671  34.09 6.7e-05   6.55
9820         TLCD2       17_2     0.891  25.23 6.5e-05  -4.79

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
9157      ZNF518B       4_11     0.000 1433.04 0.0e+00  58.12
2475       SLC2A9       4_11     0.000 1391.98 0.0e+00  53.26
2737       TRIM38       6_20     0.005  315.61 4.5e-06 -23.44
9813         MUC1       1_77     1.000  246.05 7.2e-04 -18.55
8177        THBS3       1_77     0.890  216.19 5.6e-04  16.74
187           HFE       6_20     0.005  112.31 1.7e-06 -16.15
12299   U91328.19       6_20     0.007  112.39 2.4e-06  15.64
8691      MAP3K11      11_36     0.000  166.15 3.7e-15  15.38
2953        NRBP1       2_16     0.021  215.65 1.3e-05  15.35
2956        SNX17       2_16     0.030  196.19 1.7e-05  14.39
9933       BTN3A2       6_20     0.010   88.78 2.7e-06 -13.24
12027  CTA-14H9.5       6_20     0.016   91.31 4.3e-06 -13.00
9350    HIST1H2BC       6_20     0.005   70.83 1.1e-06 -12.78
10021     ZKSCAN4       6_22     0.155  154.44 7.0e-05  12.61
11023       SIPA1      11_36     0.000  214.21 8.6e-11  12.60
6712      ZSCAN12       6_22     0.104  154.07 4.7e-05 -12.57
10765     ZDHHC18       1_18     1.000  142.25 4.1e-04  12.20
12301 RP1-86C11.7       6_21     0.000  200.97 1.7e-07 -12.06
10762      ASAH2B      10_33     0.015  128.27 5.4e-06  11.56
9293       R3HDM2      12_36     0.040  114.76 1.3e-05 -11.33

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.02631816
#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
9157      ZNF518B       4_11     0.000 1433.04 0.0e+00  58.12
2475       SLC2A9       4_11     0.000 1391.98 0.0e+00  53.26
2737       TRIM38       6_20     0.005  315.61 4.5e-06 -23.44
9813         MUC1       1_77     1.000  246.05 7.2e-04 -18.55
8177        THBS3       1_77     0.890  216.19 5.6e-04  16.74
187           HFE       6_20     0.005  112.31 1.7e-06 -16.15
12299   U91328.19       6_20     0.007  112.39 2.4e-06  15.64
8691      MAP3K11      11_36     0.000  166.15 3.7e-15  15.38
2953        NRBP1       2_16     0.021  215.65 1.3e-05  15.35
2956        SNX17       2_16     0.030  196.19 1.7e-05  14.39
9933       BTN3A2       6_20     0.010   88.78 2.7e-06 -13.24
12027  CTA-14H9.5       6_20     0.016   91.31 4.3e-06 -13.00
9350    HIST1H2BC       6_20     0.005   70.83 1.1e-06 -12.78
10021     ZKSCAN4       6_22     0.155  154.44 7.0e-05  12.61
11023       SIPA1      11_36     0.000  214.21 8.6e-11  12.60
6712      ZSCAN12       6_22     0.104  154.07 4.7e-05 -12.57
10765     ZDHHC18       1_18     1.000  142.25 4.1e-04  12.20
12301 RP1-86C11.7       6_21     0.000  200.97 1.7e-07 -12.06
10762      ASAH2B      10_33     0.015  128.27 5.4e-06  11.56
9293       R3HDM2      12_36     0.040  114.76 1.3e-05 -11.33

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: 4_11"
     genename region_tag susie_pip     mu2 PVE     z
2475   SLC2A9       4_11         0 1391.98   0 53.26
9157  ZNF518B       4_11         0 1433.04   0 58.12

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 6_20"
        genename region_tag susie_pip    mu2     PVE      z
2737      TRIM38       6_20     0.005 315.61 4.5e-06 -23.44
12299  U91328.19       6_20     0.007 112.39 2.4e-06  15.64
12332  U91328.22       6_20     0.007  13.56 2.7e-07   1.70
12483   HIST1H3A       6_20     0.005  57.42 8.6e-07  -8.42
10043   HIST1H1C       6_20     0.008  39.66 9.5e-07   4.31
187          HFE       6_20     0.005 112.31 1.7e-06 -16.15
10373   HIST1H4C       6_20     0.012  15.23 5.3e-07  -1.83
10005   HIST1H1T       6_20     0.008  84.91 1.9e-06  11.23
9350   HIST1H2BC       6_20     0.005  70.83 1.1e-06 -12.78
9348   HIST1H2AC       6_20     0.005  52.10 8.0e-07 -10.66
6686   HIST1H2BD       6_20     0.007  33.25 7.2e-07  -5.23
12425  HIST1H2BE       6_20     0.016  30.38 1.4e-06   3.69
12528  HIST1H2BF       6_20     0.005   6.80 9.8e-08  -1.88
10342  HIST1H2AD       6_20     0.011  10.83 3.4e-07  -0.12
12518   HIST1H4E       6_20     0.010  10.30 3.0e-07   0.02
12520  HIST1H2AE       6_20     0.005   6.42 1.0e-07   0.04
12444   HIST1H3E       6_20     0.021  21.39 1.3e-06   2.68
12482  HIST1H2BH       6_20     0.005  13.99 2.0e-07  -1.79
12409   HIST1H3G       6_20     0.007  19.06 3.9e-07  -5.30
6688    HIST1H4H       6_20     0.005  11.12 1.8e-07   2.02
2669      BTN3A3       6_20     0.005   6.39 9.0e-08   0.57
9933      BTN3A2       6_20     0.010  88.78 2.7e-06 -13.24
3723      BTN2A2       6_20     0.084  51.74 1.3e-05   7.22
308       BTN3A1       6_20     0.005   8.04 1.2e-07  -3.09
2771      BTN2A1       6_20     0.009  12.36 3.1e-07  -1.72
12027 CTA-14H9.5       6_20     0.016  91.31 4.3e-06 -13.00
9548       HMGN4       6_20     0.009   9.61 2.4e-07   0.93
11308      HCG11       6_20     0.008   8.72 2.0e-07  -1.60
5867        ABT1       6_20     0.056  30.31 4.9e-06   4.09
9407      ZNF322       6_20     0.007  17.98 3.5e-07   5.16

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 1_77"
          genename region_tag susie_pip    mu2     PVE      z
7201          PMVK       1_77     0.006  33.56 5.9e-07   4.16
7202        PBXIP1       1_77     0.006  33.05 5.6e-07   4.11
6909          SHC1       1_77     0.001   5.76 1.1e-08   0.78
8673         CKS1B       1_77     0.001   5.69 1.1e-08  -0.73
6907        ZBTB7B       1_77     0.002  22.16 1.3e-07  -3.16
7204         DCST2       1_77     0.058  64.79 1.1e-05   7.14
7205         DCST1       1_77     0.058  64.79 1.1e-05  -7.14
5634        ADAM15       1_77     0.001  19.76 7.6e-08  -0.59
5644         EFNA3       1_77     0.001   6.76 1.3e-08   1.16
8179         EFNA1       1_77     0.002  37.49 2.5e-07   4.74
8178       SLC50A1       1_77     0.001   7.25 1.5e-08  -1.62
8670          MTX1       1_77     0.001  81.54 1.7e-07  -9.72
9813          MUC1       1_77     1.000 246.05 7.2e-04 -18.55
8177         THBS3       1_77     0.890 216.19 5.6e-04  16.74
9103           GBA       1_77     0.001  14.72 3.6e-08   3.90
6918       FAM189B       1_77     0.002  11.67 5.6e-08  -0.20
5649          HCN3       1_77     0.001  68.05 2.9e-07   1.50
6917          FDPS       1_77     0.002  72.97 4.5e-07   0.34
4425          DAP3       1_77     0.001   8.78 1.6e-08   0.93
7207        YY1AP1       1_77     0.001  38.40 1.4e-07   5.21
4434         SYT11       1_77     0.087  90.55 2.3e-05  -1.33
5648          RIT1       1_77     0.001  25.86 5.1e-08  -1.76
4426      KIAA0907       1_77     0.002  22.14 1.1e-07  -2.57
3099       ARHGEF2       1_77     0.001  62.77 1.7e-07  -1.56
7227          SSR2       1_77     0.001  11.61 2.7e-08  -2.36
3100       LAMTOR2       1_77     0.001   7.44 1.8e-08  -0.72
4430         RAB25       1_77     0.001  23.70 4.4e-08   0.87
10224       SEMA4A       1_77     0.001  12.77 2.5e-08  -3.38
6920          PMF1       1_77     0.008  34.91 7.7e-07   3.77
11586        BGLAP       1_77     0.009  28.86 7.4e-07  -2.38
6919         PAQR6       1_77     0.003  19.18 1.7e-07  -1.66
7226        TMEM79       1_77     0.003  20.70 1.9e-07  -2.09
10714         SMG5       1_77     0.003  20.03 1.7e-07  -2.06
10641         GLMP       1_77     0.003  19.54 1.6e-07  -2.03
7224         TSACC       1_77     0.001   5.49 1.1e-08   0.50
7225          CCT3       1_77     0.001   4.92 9.4e-09  -0.21
3101         MEF2D       1_77     0.001   8.93 2.4e-08  -0.77
9652        IQGAP3       1_77     0.001   5.08 9.7e-09   0.20
10047        TTC24       1_77     0.006  22.57 3.9e-07  -2.27
7210          NAXE       1_77     0.001  11.00 4.0e-08   1.23
4428          BCAN       1_77     0.001   8.55 2.5e-08  -1.03
5586        RRNAD1       1_77     0.001   6.48 1.5e-08  -0.68
5588       ISG20L2       1_77     0.001   6.48 1.5e-08   0.68
5590          HDGF       1_77     0.001   6.75 1.5e-08   0.61
10592        NTRK1       1_77     0.001   7.80 2.0e-08  -0.81
314         SH2D2A       1_77     0.001   6.20 1.3e-08   0.52
10039        PEAR1       1_77     0.001  11.90 4.8e-08  -1.26
4429      ARHGEF11       1_77     0.001   4.94 9.3e-09   0.21
12235 RP11-85G21.3       1_77     0.001   6.58 1.5e-08   0.64
5585         FCRL5       1_77     0.001   5.17 9.9e-09  -0.21
6928         FCRL3       1_77     0.001   6.32 1.4e-08   0.79
4432         FCRL2       1_77     0.002  13.06 6.1e-08  -1.52
7243         FCRL1       1_77     0.001   7.02 1.8e-08  -1.02

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 11_36"
            genename region_tag susie_pip    mu2     PVE     z
11809          EEF1G      11_36     0.000  23.73 3.7e-15 -2.04
6091            EML3      11_36     0.000  30.93 1.0e-14 -2.40
9762           INTS5      11_36     0.000  23.89 3.8e-15 -2.05
6097          B3GAT3      11_36     0.000  21.41 2.6e-15 -1.91
1234           GANAB      11_36     0.000  73.57 2.1e-12 -3.89
11114        METTL12      11_36     0.000  33.31 1.4e-14 -2.50
12568       C11orf98      11_36     0.000  33.31 1.4e-14  2.50
10878          UQCC3      11_36     0.000  88.37 1.2e-11  4.29
7034           LBHD1      11_36     0.000  42.50 4.7e-14 -2.88
7033            GNG3      11_36     0.000   6.87 1.8e-16  0.68
8012           BSCL2      11_36     0.000  84.73 7.9e-12  4.19
7035           TTC9C      11_36     0.000  86.28 9.4e-12 -4.23
9809        TMEM179B      11_36     0.000  19.93 2.1e-15  1.83
8091         TMEM223      11_36     0.000   4.87 1.0e-16 -0.14
7037            STX5      11_36     0.000  15.30 9.8e-16  1.52
12607 RP11-727F15.14      11_36     0.000  15.25 9.7e-16  1.52
8013         HRASLS5      11_36     0.000  82.85 6.3e-12  4.14
4490         RARRES3      11_36     0.000   4.77 1.0e-16  0.03
9004         PLA2G16      11_36     0.000  33.71 1.5e-14 -2.52
9726            ATL3      11_36     0.000  11.97 5.4e-16 -1.26
4489            RTN3      11_36     0.000  14.15 8.0e-16  1.44
8014        C11orf84      11_36     0.000   7.59 2.2e-16 -0.79
803            MARK2      11_36     0.000  57.37 3.0e-13  3.40
2552           NAA40      11_36     0.000 114.95 2.5e-10  4.92
3901           PRDX5      11_36     0.000   6.77 1.8e-16  0.41
8986           COX8A      11_36     0.000  88.61 1.2e-11 -4.29
7973           OTUB1      11_36     0.000  88.68 1.2e-11 -4.29
3904           FLRT1      11_36     0.000   6.36 1.4e-16  2.37
4486         MACROD1      11_36     0.000  66.93 7.7e-13  5.16
6113           TRPT1      11_36     0.000   6.82 1.6e-16  1.40
8709           VEGFB      11_36     0.000  37.90 6.2e-15  2.93
2502          DNAJC4      11_36     0.000  14.28 5.2e-16  3.66
8708           FKBP2      11_36     0.000  81.64 3.3e-12  5.84
6114           PLCB3      11_36     0.000  38.19 3.1e-14 -4.91
15               BAD      11_36     0.000   6.89 1.7e-16  1.19
8666           ESRRA      11_36     0.000  24.31 4.4e-15  2.72
8660         TRMT112      11_36     0.000  32.42 5.8e-15 -3.51
8027         CCDC88B      11_36     0.000  22.67 7.5e-15  3.68
7041         RPS6KA4      11_36     0.000  26.66 2.2e-15 -5.96
11371     AP003774.6      11_36     0.000   5.92 1.3e-16 -2.02
9452      AP003774.4      11_36     0.000  11.56 3.8e-16 -1.23
8025             SF1      11_36     0.000   8.62 2.5e-16  1.32
4531            MEN1      11_36     0.000  65.18 9.2e-12  0.63
8432        CDC42BPG      11_36     0.000 125.35 1.1e-10 -8.15
2508            EHD1      11_36     0.000 116.13 2.5e-11 -6.82
2507           ATG2A      11_36     0.000 105.14 9.6e-12 -6.50
718          PPP2R5B      11_36     0.000 162.62 2.2e-09  8.17
8026           MAJIN      11_36     0.000 184.60 1.5e-08 -8.68
11026           ARL2      11_36     0.000 127.27 1.1e-09  4.34
8024           BATF2      11_36     0.000   9.13 2.0e-16  3.95
2504           SNX15      11_36     0.007 207.40 3.9e-06 -6.51
6116          TM7SF2      11_36     0.000  38.85 3.0e-15 -6.59
8786          ZNHIT2      11_36     0.000  38.54 2.9e-15 -6.57
6115             FAU      11_36     0.000 169.49 7.0e-08  6.98
10864          SPDYC      11_36     0.001 230.89 9.2e-07 -9.34
240            CAPN1      11_36     0.000  48.33 5.1e-14  4.29
238            POLA2      11_36     0.000   9.30 2.2e-16  3.75
4530            DPF2      11_36     0.000  10.75 5.4e-16  2.49
8745           TIGD3      11_36     0.000  17.13 1.5e-15 -1.57
7038        SLC25A45      11_36     0.000  10.39 5.2e-16 -1.80
11663          NEAT1      11_36     0.000  29.53 2.5e-14 -1.22
8022           LTBP3      11_36     0.000  18.39 3.1e-15 -1.32
9040          FAM89B      11_36     0.000  11.21 9.9e-16 -0.29
8693           KCNK7      11_36     0.000 114.27 7.4e-14  9.88
8691         MAP3K11      11_36     0.000 166.15 3.7e-15 15.38
10388          PCNX3      11_36     0.000  45.87 7.7e-14  5.06
11023          SIPA1      11_36     0.000 214.21 8.6e-11 12.60
8652            RELA      11_36     0.000  47.96 3.0e-15  5.28
8647            KAT5      11_36     0.000  65.35 2.1e-12 -4.87
8636        RNASEH2C      11_36     0.000  77.38 1.7e-12  6.70
11839   RP11-770G2.2      11_36     0.000  37.20 8.3e-16  6.82
11803          AP5B1      11_36     0.000  24.81 1.1e-14  0.42
8610           MUS81      11_36     0.000  63.38 3.0e-12 -0.55
8613            CFL1      11_36     0.000  57.18 9.0e-13 -0.96
8599          EFEMP2      11_36     0.000   7.74 1.6e-16 -1.54
8590            CTSW      11_36     0.000  22.38 1.4e-15 -2.31
8586            FIBP      11_36     0.000  85.14 8.9e-11  1.45
8924         CCDC85B      11_36     0.000  22.56 6.5e-16  4.98
8918        C11orf68      11_36     0.000  44.46 8.1e-15  5.40
8915           DRAP1      11_36     0.000  55.80 2.4e-15  8.69
8909        TSGA10IP      11_36     0.000  33.16 1.0e-15  7.26
8890           BANF1      11_36     0.000  60.39 1.3e-15  9.54
8888            CST6      11_36     0.000  24.19 4.8e-15  1.02
1166           SF3B2      11_36     0.000  16.65 4.7e-16  3.24
8863           PACS1      11_36     0.000  41.86 6.2e-15 -5.26

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 2_16"
      genename region_tag susie_pip    mu2     PVE     z
11045  SLC35F6       2_16     0.015  18.26 7.8e-07 -3.70
3366   TMEM214       2_16     0.015   7.99 3.5e-07  1.10
5074   EMILIN1       2_16     0.021  29.91 1.8e-06  4.84
5061       KHK       2_16     0.013   6.76 2.6e-07 -0.99
5059    CGREF1       2_16     0.011   4.82 1.5e-07 -0.06
5070      PREB       2_16     0.011   8.91 3.0e-07 -2.00
5076    ATRAID       2_16     0.018  77.61 4.0e-06 -9.11
1090       CAD       2_16     0.046  21.48 2.9e-06 -3.07
5071    SLC5A6       2_16     0.018  16.85 8.9e-07  2.99
7303       UCN       2_16     0.016  20.68 9.6e-07 -4.41
2952    GTF3C2       2_16     0.017  21.62 1.1e-06  4.48
2956     SNX17       2_16     0.030 196.19 1.7e-05 14.39
7304    ZNF513       2_16     0.011  96.76 3.2e-06 10.03
2953     NRBP1       2_16     0.021 215.65 1.3e-05 15.35
5057    IFT172       2_16     0.027  32.83 2.5e-06 -5.41
1087      GCKR       2_16     0.029  34.30 2.9e-06  5.50
10613     GPN1       2_16     0.029  43.15 3.6e-06  6.44
9018   CCDC121       2_16     0.187  27.27 1.5e-05  1.38
6660       BRE       2_16     0.038  30.25 3.3e-06  4.67

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
56720  rs766167074      1_118     1.000 5751.16 1.7e-02  -2.94
72399     rs780093       2_16     1.000  441.15 1.3e-03 -22.26
100027    rs882759       2_70     1.000   56.04 1.6e-04  -7.38
113237   rs7565788      2_103     1.000  140.49 4.1e-04 -11.91
113299   rs6433115      2_103     1.000   50.38 1.5e-04   7.40
115592    rs863678      2_106     1.000  102.08 3.0e-04   6.01
117533  rs11690832      2_110     1.000   59.14 1.7e-04   8.04
126605   rs2068218      2_128     1.000   44.83 1.3e-04  -6.60
198856 rs141435299       4_10     1.000  514.28 1.5e-03  -1.03
199060  rs57136958       4_11     1.000  377.71 1.1e-03 -10.08
199116  rs13115469       4_11     1.000 9484.50 2.8e-02 133.38
199151   rs3775948       4_11     1.000 9730.85 2.8e-02 131.05
199177  rs75968456       4_11     1.000  660.26 1.9e-03  -2.69
199480   rs6831973       4_12     1.000   85.49 2.5e-04  -9.83
199506  rs76285604       4_12     1.000  134.17 3.9e-04 -11.85
199593 rs142309009       4_12     1.000   73.37 2.1e-04  -9.06
225352 rs149027545       4_59     1.000 2233.20 6.5e-03  53.88
225418  rs10022462       4_60     1.000   99.49 2.9e-04 -12.99
229800  rs35518360       4_67     1.000   45.37 1.3e-04  -6.68
276354    rs255749       5_31     1.000   71.94 2.1e-04   8.64
277125  rs10077826       5_33     1.000   35.87 1.0e-04  -5.73
281562  rs10942549       5_43     1.000  170.65 5.0e-04 -15.21
318371    rs630258        6_7     1.000  163.02 4.7e-04 -17.15
324838   rs7754961       6_19     1.000  123.00 3.6e-04  15.01
324877  rs12213398       6_19     1.000   45.86 1.3e-04  -4.69
325521   rs6908155       6_21     1.000   63.06 1.8e-04  -2.88
329970   rs2856992       6_27     1.000   50.60 1.5e-04  -5.75
358655  rs10782229       6_84     1.000   42.78 1.2e-04   5.77
373878  rs78148157        7_3     1.000  114.56 3.3e-04  -9.55
373879  rs13241427        7_3     1.000   86.95 2.5e-04   8.10
400322 rs761767938       7_49     1.000 6271.03 1.8e-02  -3.88
400330   rs1544459       7_49     1.000 6306.74 1.8e-02  -4.00
449033   rs2941484       8_56     1.000  159.21 4.6e-04  15.49
483652  rs56030777       9_25     1.000  105.17 3.1e-04  -4.38
497215   rs1800977       9_53     1.000   48.11 1.4e-04  -6.96
524162  rs35182775      10_33     1.000  204.88 6.0e-04  15.09
527297  rs11510917      10_39     1.000  214.66 6.2e-04 -19.07
527310   rs1171619      10_39     1.000  332.94 9.7e-04  21.17
534222 rs149429992      10_52     1.000 7903.90 2.3e-02   2.53
545268  rs10886117      10_72     1.000   54.68 1.6e-04   7.40
562711 rs369062552      11_21     1.000  184.91 5.4e-04  11.72
562721  rs34830202      11_21     1.000  193.74 5.6e-04 -12.09
598886  rs11056397      12_13     1.000   35.29 1.0e-04  -5.84
611235   rs7397189      12_36     1.000  338.53 9.8e-04 -20.09
627821    rs653178      12_67     1.000  169.46 4.9e-04 -14.16
630757   rs2701194      12_74     1.000   71.82 2.1e-04   8.03
635090  rs76734539      12_82     1.000   60.40 1.8e-04   7.62
647865 rs566812111      13_25     1.000 2967.96 8.6e-03   2.56
677618  rs72681869      14_20     1.000   65.36 1.9e-04  -8.26
711310   rs2472297      15_35     1.000   66.83 1.9e-04  -9.72
711552 rs145727191      15_35     1.000   51.54 1.5e-04   8.91
711581   rs2955742      15_36     1.000   42.84 1.2e-04   7.22
719058  rs59646751      15_48     1.000  149.60 4.4e-04  14.47
732653  rs12927956      16_27     1.000   69.00 2.0e-04   7.55
760194   rs3794748      17_32     1.000  171.13 5.0e-04 -13.59
802776  rs56287732      19_23     1.000   40.84 1.2e-04   5.41
824742    rs209955      20_32     1.000   36.46 1.1e-04   5.89
836541    rs219783      21_16     1.000   88.59 2.6e-04  -9.60
868178  rs75246752       1_73     1.000  308.73 9.0e-04 -18.17
868325 rs185073199       1_73     1.000  281.42 8.2e-04 -18.84
930928 rs140927145       7_92     1.000 8892.68 2.6e-02  -2.70
940154 rs145700013       7_94     1.000   58.23 1.7e-04  -2.96
957040 rs542984928      11_36     1.000  241.27 7.0e-04  23.70
969052   rs6581124      12_35     1.000  109.18 3.2e-04 -10.40
12937   rs10788884       1_30     0.999   41.08 1.2e-04   6.38
507175  rs12380852       9_73     0.999   34.36 1.0e-04   5.71
552287   rs3842748       11_2     0.999   85.12 2.5e-04  -8.12
647869  rs12430288      13_25     0.999 2995.88 8.7e-03   2.63
755421      rs2688      17_22     0.999   32.30 9.4e-05  -5.49
957069  rs12363578      11_36     0.999  495.30 1.4e-03 -26.87
199276   rs4140694       4_11     0.998  864.91 2.5e-03  16.31
199599  rs61795273       4_12     0.998   54.59 1.6e-04  -7.25
534224   rs2152629      10_52     0.998 7922.61 2.3e-02   2.47
571707   rs1203614      11_37     0.998   50.15 1.5e-04   6.02
761515  rs11650989      17_36     0.998   39.50 1.1e-04   7.88
775289    rs527616      18_14     0.998   31.04 9.0e-05  -5.57
806242    rs814573      19_31     0.998   30.49 8.8e-05  -5.56
243285   rs4521364       4_95     0.997   56.37 1.6e-04  -5.84
314173 rs139078584      5_106     0.997   31.64 9.2e-05   6.53
617714   rs1848968      12_48     0.997   45.96 1.3e-04  -6.75
318194 rs200823080        6_6     0.996   33.74 9.8e-05   5.68
551072   rs2767419      10_85     0.996   34.60 1.0e-04  -5.65
179427 rs145422957       3_92     0.995   31.62 9.2e-05  -5.47
763623 rs113408695      17_39     0.995   31.74 9.2e-05  -5.55
204001    rs358256       4_20     0.994   33.89 9.8e-05   5.69
325354  rs13191326       6_21     0.994  231.67 6.7e-04  13.59
930924   rs6954405       7_92     0.994 8872.27 2.6e-02  -2.37
86985    rs7561263       2_46     0.993   35.22 1.0e-04  -6.03
173059  rs11712964       3_78     0.993   30.33 8.8e-05   5.45
713099   rs7174325      15_38     0.993   28.69 8.3e-05   4.89
826463   rs1407040      20_34     0.993   29.66 8.6e-05  -5.21
84547     rs778147       2_40     0.992   48.03 1.4e-04   6.79
512387  rs72777711      10_10     0.991   30.41 8.8e-05  -5.28
802630  rs71176165      19_23     0.991   69.48 2.0e-04  -9.09
719065   rs8037855      15_48     0.990   64.88 1.9e-04  10.82
92136   rs10196697       2_54     0.989   32.94 9.5e-05  -5.64
406723  rs45467892       7_61     0.989   41.19 1.2e-04  -6.54
715842 rs113404146      15_43     0.989   34.94 1.0e-04   5.81
400326  rs11972122       7_49     0.988 5827.76 1.7e-02  -3.92
692444  rs73349296      14_50     0.988   43.25 1.2e-04  -6.59
449071   rs2941432       8_56     0.987   69.94 2.0e-04  11.61
177014 rs115604285       3_87     0.986   32.60 9.3e-05   6.16
684650  rs10151620      14_34     0.985   36.60 1.0e-04   6.06
589995  rs73022311      11_77     0.984   26.57 7.6e-05  -4.89
738367  rs72799826      16_38     0.984   32.73 9.4e-05  -5.84
391262    rs700752       7_34     0.983   48.42 1.4e-04   6.67
844931    rs740219      22_10     0.982   35.74 1.0e-04   6.06
230575   rs2903386       4_69     0.979   34.27 9.8e-05   5.82
540386   rs7900763      10_64     0.978   26.87 7.6e-05  -4.72
394783 rs140420256       7_39     0.976   25.43 7.2e-05   4.74
609668   rs1878234      12_31     0.976   29.89 8.5e-05  -5.23
210003  rs12639940       4_32     0.974   25.93 7.3e-05  -4.41
330172   rs1126511       6_27     0.974   50.54 1.4e-04   5.39
630766  rs80019595      12_74     0.974   54.21 1.5e-04  -6.72
802814    rs889140      19_23     0.971   29.70 8.4e-05  -4.20
572155   rs6591334      11_37     0.970   38.19 1.1e-04  -5.74
712955   rs8024096      15_37     0.970   28.38 8.0e-05   5.11
992645   rs3760994       19_2     0.970   46.54 1.3e-04  -6.66
150789  rs62259692       3_36     0.969   35.30 9.9e-05   7.65
204798  rs34811474       4_21     0.969   35.24 9.9e-05  -5.82
429068  rs17151140       8_13     0.969   55.65 1.6e-04   5.11
738028  rs56259873      16_37     0.969   90.29 2.5e-04  10.62
775330    rs162000      18_14     0.968   25.57 7.2e-05   4.95
237565  rs62323480       4_83     0.966   25.77 7.2e-05   4.40
247282 rs139212650      4_102     0.965   24.93 7.0e-05   4.61
6284    rs10917555       1_13     0.964   25.92 7.3e-05  -4.75
34463    rs1979575       1_75     0.964   25.34 7.1e-05  -4.39
422135  rs11761498       7_98     0.961   32.78 9.2e-05   5.26
479025  rs10964603       9_16     0.961   25.17 7.0e-05  -4.52
277286 rs536916238       5_33     0.959   31.98 8.9e-05  -0.30
922914  rs56142449       6_61     0.958   47.74 1.3e-04   7.05
438509 rs111965375       8_34     0.955   25.87 7.2e-05   5.22
570268 rs199725747      11_32     0.955   31.05 8.6e-05   5.34
225372  rs28366540       4_59     0.952  558.20 1.5e-03 -33.87
333088  rs10223666       6_34     0.952  191.43 5.3e-04  14.34
70389   rs62112223       2_12     0.951   30.04 8.3e-05  -5.32
826332  rs73185042      20_34     0.950   25.34 7.0e-05  -4.69
316516   rs1272694        6_3     0.949   34.82 9.6e-05  -5.95
798130  rs11668601      19_14     0.948   34.35 9.5e-05   7.11
326301  rs13437375       6_23     0.944   29.96 8.2e-05   2.72
49268  rs114934802      1_104     0.941   23.80 6.5e-05  -3.09
277268    rs173964       5_33     0.941  104.26 2.9e-04   8.09
463617    rs921719       8_83     0.940   27.15 7.4e-05  -5.17
868345 rs139428292       1_73     0.940   75.54 2.1e-04  -9.43
140620   rs6778028       3_12     0.937   25.47 6.9e-05   4.59
424506 rs117950418        8_4     0.937   24.64 6.7e-05  -4.62
683084  rs34528648      14_33     0.935   32.20 8.8e-05   5.45
721644   rs2601781       16_4     0.933   25.40 6.9e-05   4.71
719041  rs75422555      15_47     0.932   31.13 8.4e-05  -6.41
439799  rs12543287       8_37     0.931   34.57 9.4e-05  -5.67
157749  rs56324130       3_49     0.927   23.89 6.4e-05  -4.43
218953 rs144193499       4_48     0.927   25.06 6.8e-05   4.55
505329 rs201421930       9_69     0.924   34.47 9.3e-05  -5.88
551018  rs75184896      10_84     0.922   26.72 7.2e-05   5.47
193668   rs7642977      3_118     0.921   29.16 7.8e-05   5.14
352883  rs12196331       6_71     0.921   28.67 7.7e-05   5.24
718742 rs116887089      15_47     0.920   24.30 6.5e-05   4.45
236735  rs72680231       4_81     0.919   24.59 6.6e-05  -4.66
185840   rs1290790      3_104     0.916   69.97 1.9e-04   6.15
610857  rs11608918      12_33     0.916   26.13 7.0e-05   4.71
235872  rs41278087       4_79     0.915   24.02 6.4e-05  -4.60
464936  rs36041912       8_85     0.911   24.43 6.5e-05  -4.56
522457  rs58434594      10_30     0.911   27.86 7.4e-05   4.83
154979  rs56145049       3_43     0.910   24.81 6.6e-05  -4.67
271612  rs13170671       5_23     0.909   70.56 1.9e-04   8.54
315427   rs6873880      5_109     0.904   27.94 7.3e-05  -5.03
346960 rs118165878       6_58     0.901   23.85 6.2e-05   4.44
211738 rs112161979       4_35     0.899   23.68 6.2e-05   4.56
229866  rs13140033       4_68     0.898   24.39 6.4e-05  -4.46
66812  rs139638572        2_6     0.896   28.26 7.4e-05   4.44
150972   rs2276816       3_36     0.892   32.58 8.5e-05  -0.96
72596    rs7606480       2_19     0.891   26.38 6.8e-05  -4.89
737597  rs12599580      16_36     0.891   27.48 7.1e-05   5.66
783917    rs784218      18_30     0.891   27.10 7.0e-05   4.29
711522 rs553274247      15_35     0.890  122.43 3.2e-04  -5.27
280971   rs6886422       5_42     0.886   25.31 6.5e-05  -4.52
447667  rs34650318       8_53     0.884   24.60 6.3e-05  -4.51
189108  rs17593458      3_110     0.880   25.58 6.5e-05   4.52
287010   rs3952745       5_53     0.874   32.71 8.3e-05  -6.44
8429     rs7516039       1_20     0.873   25.70 6.5e-05   4.71
324791  rs62392365       6_19     0.873   72.25 1.8e-04  11.55
910347  rs78595810       3_33     0.873   37.99 9.6e-05   6.24
354984   rs1997649       6_76     0.870   23.80 6.0e-05   4.28
572496  rs10752584      11_38     0.867   30.20 7.6e-05   4.73
613671 rs113897089      12_40     0.866   29.67 7.5e-05   5.22
544001  rs11196217      10_70     0.864   25.06 6.3e-05  -4.56
198439  rs28680668        4_9     0.863   33.18 8.3e-05   5.66
810362   rs2003700       20_1     0.863   24.39 6.1e-05   4.41
552274  rs11042594       11_2     0.862   67.85 1.7e-04   6.63
956689  rs78915580      11_36     0.862   65.08 1.6e-04  -7.09
792984   rs2074457       19_4     0.861   25.90 6.5e-05   4.44
5089     rs4336844       1_11     0.857   34.05 8.5e-05   4.45
571730    rs509533      11_37     0.855   45.23 1.1e-04   7.93
959380   rs5792371      11_36     0.853  267.54 6.6e-04  18.90
413023  rs75520353       7_74     0.852   24.12 6.0e-05   4.36
227402  rs11722010       4_63     0.851   25.08 6.2e-05  -4.49
698690  rs62004133       15_8     0.851   25.21 6.2e-05  -4.56
741707  rs76862947      16_44     0.849  139.45 3.4e-04 -11.76
26425    rs9432440       1_58     0.845   25.44 6.3e-05   4.55
483621  rs10971408       9_25     0.845   28.05 6.9e-05  -3.63
99801    rs2311597       2_70     0.838   59.90 1.5e-04  -9.93
40871     rs604388       1_87     0.836   24.72 6.0e-05   4.43
446348  rs71553284       8_50     0.836   25.43 6.2e-05   4.39
785487  rs17773471      18_33     0.835   35.77 8.7e-05   7.37
103829  rs71862935       2_79     0.832   27.40 6.6e-05  -4.89
737963  rs72799382      16_37     0.832   42.87 1.0e-04  -6.78
110874  rs12467636       2_96     0.830   35.47 8.6e-05   5.79
122243  rs56059523      2_120     0.829   29.74 7.2e-05   4.32
302294    rs356486       5_82     0.828   29.34 7.1e-05   5.01
150575  rs78342753       3_35     0.824   34.63 8.3e-05   5.19
407387   rs7803747       7_63     0.824   26.62 6.4e-05   5.00
601668   rs4077798      12_18     0.824   30.58 7.3e-05  -5.30
399594   rs1179610       7_48     0.820   26.15 6.2e-05   4.72
805138  rs11671669      19_30     0.816   26.29 6.2e-05   4.49
580024  rs57569860      11_52     0.815   23.94 5.7e-05   3.53
139938   rs9846767       3_11     0.811   25.13 5.9e-05   4.59
87000   rs72837990       2_46     0.810   24.87 5.9e-05  -4.42
435464   rs4871905       8_24     0.810  131.51 3.1e-04  12.10
824772 rs570322378      20_32     0.809   27.77 6.5e-05   4.45
57987     rs678615      1_122     0.807   40.29 9.5e-05  -6.17
640192   rs9315009       13_9     0.807   36.78 8.6e-05   5.93
547754  rs10794177      10_78     0.806   52.00 1.2e-04   7.11
547877   rs1693628      10_78     0.802   26.48 6.2e-05  -4.38

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
199151   rs3775948       4_11     1.000 9730.85 2.8e-02  131.05
199118   rs6838021       4_11     0.000 9648.29 0.0e+00  133.98
199119   rs6823324       4_11     0.000 9629.18 0.0e+00 -133.15
199122  rs11723439       4_11     0.000 9489.36 0.0e+00 -132.41
199116  rs13115469       4_11     1.000 9484.50 2.8e-02  133.38
930928 rs140927145       7_92     1.000 8892.68 2.6e-02   -2.70
930924   rs6954405       7_92     0.994 8872.27 2.6e-02   -2.37
930923   rs6953180       7_92     0.299 8870.22 7.7e-03   -2.31
930922  rs11971515       7_92     0.090 8868.53 2.3e-03   -2.27
930921  rs11974627       7_92     0.101 8867.27 2.6e-03   -2.27
930926 rs528981137       7_92     0.006 8864.62 1.6e-04   -2.16
930927 rs549027364       7_92     0.007 8864.41 1.7e-04   -2.16
199141   rs7439210       4_11     0.000 8719.19 0.0e+00 -128.12
930918   rs7807788       7_92     0.000 8649.45 8.7e-08   -2.43
930912  rs11546289       7_92     0.000 8625.48 1.9e-09   -2.37
930913  rs13243678       7_92     0.000 8616.95 7.4e-10   -2.36
930914  rs11546290       7_92     0.000 8611.36 2.4e-10   -2.34
930910   rs6952916       7_92     0.000 8570.29 4.4e-12   -2.29
930907   rs7793916       7_92     0.000 8558.99 3.0e-12   -2.30
930908   rs7781312       7_92     0.000 8539.69 4.1e-13   -2.28
930906  rs13247593       7_92     0.000 8530.79 8.3e-12   -2.38
930909   rs6965276       7_92     0.000 8484.08 6.6e-15   -2.28
930931 rs112579924       7_92     0.000 8444.92 0.0e+00   -1.71
930905 rs568222600       7_92     0.000 8406.53 3.5e-17   -2.28
930936  rs12667507       7_92     0.000 8396.38 0.0e+00    2.10
930903 rs112799047       7_92     0.000 8281.08 0.0e+00   -2.31
930904  rs10233906       7_92     0.000 8270.55 0.0e+00   -2.31
930901 rs112684174       7_92     0.000 8261.81 0.0e+00   -2.32
930902 rs112634954       7_92     0.000 8253.40 0.0e+00   -2.32
930896   rs6965440       7_92     0.000 8240.31 0.0e+00   -2.33
930895   rs6464930       7_92     0.000 8228.70 0.0e+00   -2.35
930890   rs6978068       7_92     0.000 8187.98 0.0e+00   -2.28
930881  rs10233535       7_92     0.000 8087.97 0.0e+00   -2.31
930882  rs56932055       7_92     0.000 8085.60 0.0e+00   -2.32
930886   rs6958855       7_92     0.000 8084.60 0.0e+00   -2.34
930889   rs6963547       7_92     0.000 8081.91 0.0e+00   -2.37
930871  rs10269104       7_92     0.000 8080.73 0.0e+00   -2.30
930887   rs1551926       7_92     0.000 8058.67 0.0e+00    2.45
199209  rs11723742       4_11     0.000 7995.98 0.0e+00 -116.64
930915   rs2306169       7_92     0.000 7930.34 0.0e+00    2.54
534224   rs2152629      10_52     0.998 7922.61 2.3e-02    2.47
930917   rs2306170       7_92     0.000 7913.74 0.0e+00    2.33
534220   rs7913261      10_52     0.562 7905.82 1.3e-02    2.51
534222 rs149429992      10_52     1.000 7903.90 2.3e-02    2.53
199293  rs17389602       4_11     0.000 7847.08 0.0e+00 -117.34
199295  rs78917351       4_11     0.000 7835.51 0.0e+00 -117.31
534225   rs1360953      10_52     0.000 7825.99 1.6e-08    2.50
534230  rs76574695      10_52     0.000 7824.36 1.4e-07    2.40
930925   rs6953222       7_92     0.000 7824.06 0.0e+00   -1.43
930920  rs34909003       7_92     0.000 7772.17 0.0e+00   -1.86

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
199116  rs13115469       4_11     1.000 9484.50 0.02800 133.38
199151   rs3775948       4_11     1.000 9730.85 0.02800 131.05
930924   rs6954405       7_92     0.994 8872.27 0.02600  -2.37
930928 rs140927145       7_92     1.000 8892.68 0.02600  -2.70
534222 rs149429992      10_52     1.000 7903.90 0.02300   2.53
534224   rs2152629      10_52     0.998 7922.61 0.02300   2.47
400322 rs761767938       7_49     1.000 6271.03 0.01800  -3.88
400330   rs1544459       7_49     1.000 6306.74 0.01800  -4.00
56720  rs766167074      1_118     1.000 5751.16 0.01700  -2.94
400326  rs11972122       7_49     0.988 5827.76 0.01700  -3.92
534220   rs7913261      10_52     0.562 7905.82 0.01300   2.51
647869  rs12430288      13_25     0.999 2995.88 0.00870   2.63
647865 rs566812111      13_25     1.000 2967.96 0.00860   2.56
930923   rs6953180       7_92     0.299 8870.22 0.00770  -2.31
225352 rs149027545       4_59     1.000 2233.20 0.00650  53.88
56718    rs2486737      1_118     0.370 5713.07 0.00620  -3.17
56719     rs971534      1_118     0.358 5713.05 0.00590  -3.16
56717   rs10489611      1_118     0.304 5712.90 0.00500  -3.16
56714    rs2790891      1_118     0.290 5712.52 0.00480  -3.16
56715    rs2491405      1_118     0.290 5712.52 0.00480  -3.16
56711    rs2256908      1_118     0.277 5712.49 0.00460  -3.16
647868   rs1579715      13_25     0.388 2965.08 0.00330  -2.77
56726    rs2248646      1_118     0.181 5710.15 0.00300  -3.14
56727    rs2211176      1_118     0.173 5710.35 0.00290  -3.14
56728    rs2790882      1_118     0.173 5710.35 0.00290  -3.14
930921  rs11974627       7_92     0.101 8867.27 0.00260  -2.27
199276   rs4140694       4_11     0.998  864.91 0.00250  16.31
930922  rs11971515       7_92     0.090 8868.53 0.00230  -2.27
199177  rs75968456       4_11     1.000  660.26 0.00190  -2.69
56707    rs1076804      1_118     0.089 5703.74 0.00150  -3.14
198856 rs141435299       4_10     1.000  514.28 0.00150  -1.03
225372  rs28366540       4_59     0.952  558.20 0.00150 -33.87
957069  rs12363578      11_36     0.999  495.30 0.00140 -26.87
72399     rs780093       2_16     1.000  441.15 0.00130 -22.26
56729    rs1416913      1_118     0.071 5703.06 0.00120  -3.13
56723    rs2739509      1_118     0.070 5609.25 0.00110  -3.30
199060  rs57136958       4_11     1.000  377.71 0.00110 -10.08
198851 rs146530806       4_10     0.572  606.87 0.00100   6.59
198849 rs144362537       4_10     0.558  607.49 0.00099   6.57
611235   rs7397189      12_36     1.000  338.53 0.00098 -20.09
527310   rs1171619      10_39     1.000  332.94 0.00097  21.17
198850  rs34658640       4_10     0.536  607.42 0.00095   6.57
868178  rs75246752       1_73     1.000  308.73 0.00090 -18.17
56732    rs2790874      1_118     0.051 5702.25 0.00085  -3.12
868325 rs185073199       1_73     1.000  281.42 0.00082 -18.84
957040 rs542984928      11_36     1.000  241.27 0.00070  23.70
325354  rs13191326       6_21     0.994  231.67 0.00067  13.59
959380   rs5792371      11_36     0.853  267.54 0.00066  18.90
527297  rs11510917      10_39     1.000  214.66 0.00062 -19.07
524162  rs35182775      10_33     1.000  204.88 0.00060  15.09

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
199118   rs6838021       4_11         0 9648.29 0.000  133.98
199116  rs13115469       4_11         1 9484.50 0.028  133.38
199119   rs6823324       4_11         0 9629.18 0.000 -133.15
199122  rs11723439       4_11         0 9489.36 0.000 -132.41
199151   rs3775948       4_11         1 9730.85 0.028  131.05
199141   rs7439210       4_11         0 8719.19 0.000 -128.12
199293  rs17389602       4_11         0 7847.08 0.000 -117.34
199295  rs78917351       4_11         0 7835.51 0.000 -117.31
199209  rs11723742       4_11         0 7995.98 0.000 -116.64
199347  rs11722185       4_11         0 7644.92 0.000 -115.64
199333  rs11727390       4_11         0 7650.25 0.000 -115.57
199339   rs4697745       4_11         0 7032.95 0.000 -110.85
199370 rs546391476       4_11         0 6948.16 0.000 -110.51
199321  rs10489070       4_11         0 6827.95 0.000 -109.37
199159   rs9291642       4_11         0 5586.71 0.000  103.74
199256   rs4697717       4_11         0 6463.32 0.000 -103.72
199278    rs887732       4_11         0 6247.19 0.000 -103.53
199181   rs7349721       4_11         0 5360.29 0.000  100.68
199125   rs7375643       4_11         0 4136.86 0.000  -90.99
199145  rs11723970       4_11         0 4130.56 0.000   90.75
199124   rs7375599       4_11         0 4108.12 0.000  -90.71
199377   rs5856057       4_11         0 4544.79 0.000  -86.60
199139   rs6449177       4_11         0 3551.11 0.000  -85.49
199152  rs34501273       4_11         0 3886.64 0.000   84.99
199154   rs3733586       4_11         0 3887.21 0.000   84.99
199155  rs35438220       4_11         0 3886.06 0.000   84.99
199158  rs12507330       4_11         0 3885.14 0.000   84.97
199149  rs17187075       4_11         0 3737.74 0.000   84.09
199121  rs12498956       4_11         0 3434.53 0.000  -83.90
199163   rs3756236       4_11         0 3719.29 0.000   83.79
199178  rs11727199       4_11         0 3719.96 0.000   83.20
199179  rs10939665       4_11         0 3717.49 0.000   83.17
199180  rs12508991       4_11         0 3718.21 0.000   83.17
199192  rs35501905       4_11         0 3715.15 0.000   83.12
199113  rs35955619       4_11         0 3189.77 0.000  -82.89
199131   rs6844316       4_11         0 3332.72 0.000  -82.84
199136   rs4295261       4_11         0 3334.21 0.000  -82.84
199129  rs34297373       4_11         0 3331.51 0.000  -82.82
199132  rs28837683       4_11         0 3327.14 0.000  -82.81
199138   rs7434391       4_11         0 3334.52 0.000  -82.79
199175  rs13148371       4_11         0 3696.12 0.000   82.73
199123   rs7376155       4_11         0 3308.95 0.000  -82.67
199127   rs4314284       4_11         0 3308.99 0.000  -82.67
199196   rs7678211       4_11         0 3824.66 0.000   82.60
199217   rs4235356       4_11         0 3341.83 0.000  -80.91
199166   rs6827785       4_11         0 3247.91 0.000   79.26
199191  rs13120348       4_11         0 3236.71 0.000   78.51
199183  rs13122689       4_11         0 3227.14 0.000   78.50
199184  rs12504565       4_11         0 3227.57 0.000   78.50
199176  rs12506122       4_11         0 3240.61 0.000   78.48

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] 15
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)
ZNF276 gene(s) from the input list not found in DisGeNET CURATEDZDHHC18 gene(s) from the input list not found in DisGeNET CURATEDANO7 gene(s) from the input list not found in DisGeNET CURATEDTLCD2 gene(s) from the input list not found in DisGeNET CURATEDKCNK17 gene(s) from the input list not found in DisGeNET CURATEDAGAP3 gene(s) from the input list not found in DisGeNET CURATEDCSNK1G2 gene(s) from the input list not found in DisGeNET CURATED
                                Description        FDR Ratio  BgRatio
61        Medullary cystic kidney disease 1 0.02061431   1/8   1/9703
64 Mental Retardation, Autosomal Dominant 1 0.02061431   1/8   1/9703
72            2q23.1 microdeletion syndrome 0.02061431   1/8   1/9703
4                       Cannabis Dependence 0.07375783   1/8  17/9703
5            Neoplastic Cell Transformation 0.07375783   2/8 139/9703
22                    Neoplasm Invasiveness 0.07375783   2/8 184/9703
25                     Peritoneal Neoplasms 0.07375783   1/8  10/9703
31                            Gastric ulcer 0.07375783   1/8  18/9703
32                   Aganglionosis, Colonic 0.07375783   1/8  15/9703
48      Carcinomatosis of peritoneal cavity 0.07375783   1/8  10/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