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-30650_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30660_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30670_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30680_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30690_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30700_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30710_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30720_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30730_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30740_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30750_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30760_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30770_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30780_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30790_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30800_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30810_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30820_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30830_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30840_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30850_irnt_Liver.Rmd
    Modified:   analysis/ukb-d-30860_irnt_Liver.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-30630_irnt_Whole_Blood.Rmd) and HTML (docs/ukb-d-30630_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 05a98b7 wesleycrouse 2021-08-07 adding additional results
html 05a98b7 wesleycrouse 2021-08-07 adding additional results
html 03e541c wesleycrouse 2021-07-29 Cleaning up report generation
Rmd 276893d wesleycrouse 2021-07-29 Updating reports
html 276893d wesleycrouse 2021-07-29 Updating reports

Overview

These are the results of a ctwas analysis of the UK Biobank trait Apoliprotein A (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-30630_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.0196334837 0.0001842734 
#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.11155 16.68317 
#report sample size
print(sample_size)
[1] 313387
#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.01745490 0.08531908 
#compare sum(PIP*mu2/sample_size) with above PVE calculation
c(sum(ctwas_gene_res$PVE),sum(ctwas_snp_res$PVE))
[1] 0.1325020 0.8951585

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    76.68 2.4e-04 -10.10
1980          FCGRT      19_34     1.000 29303.85 9.4e-02  -4.89
1652          PCIF1      20_28     0.999    64.50 2.1e-04   8.15
4564          PSRC1       1_67     0.998   231.34 7.4e-04  16.64
7654          PSMC3      11_29     0.997   109.99 3.5e-04 -14.44
4610           ACP2      11_29     0.996   105.08 3.3e-04 -14.33
12237 RP11-110I1.14      11_71     0.994    54.36 1.7e-04   8.05
8863          PACS1      11_36     0.986    28.53 9.0e-05   5.21
9863          LAMP1      13_62     0.980    25.00 7.8e-05   4.79
8482          SPSB1        1_6     0.979    23.06 7.2e-05   4.29
5389           CTRL      16_36     0.979   258.47 8.1e-04  15.96
3378           GPAM      10_70     0.973    99.00 3.1e-04  10.46
12304  RP11-54O7.17        1_1     0.965    43.82 1.3e-04  -6.74
6590          NTAN1      16_15     0.960    36.78 1.1e-04  -5.73
1420         MARCH2       19_8     0.955    27.63 8.4e-05  -7.04
5397          VPS53       17_1     0.952    23.82 7.2e-05   4.58
23             M6PR       12_9     0.947    35.97 1.1e-04   5.68
6370          CEBPG      19_23     0.946    24.27 7.3e-05  -5.66
8628            SP3      2_105     0.935    24.51 7.3e-05   2.76
5834        TNFAIP8       5_72     0.931    25.75 7.6e-05  -4.81
4316          MAP1B       5_43     0.920    21.79 6.4e-05   4.37
6439         SLFN13      17_21     0.920    20.51 6.0e-05   4.09
4360          TRIM5       11_4     0.912    26.92 7.8e-05   4.08
6089          FADS1      11_34     0.910   156.08 4.5e-04 -14.88
9777         RAB11B       19_8     0.906    46.72 1.4e-04 -10.38
2333         PITRM1       10_4     0.899    19.22 5.5e-05  -3.87
9322             F2      11_28     0.898    65.51 1.9e-04 -10.46
2891           GBE1       3_54     0.884    21.80 6.2e-05   4.39
2410            MLX      17_25     0.882    55.30 1.6e-04  -7.32
7513          FOXK1        7_6     0.878    25.40 7.1e-05  -4.80
4283          TRAF3      14_54     0.875    26.76 7.5e-05  -4.95
10417        MRPL21      11_38     0.874    52.78 1.5e-04   7.15
2442          IFT20      17_18     0.873    86.29 2.4e-04  -9.35
1267         PABPC4       1_24     0.868   185.89 5.1e-04  14.60
2388           BLMH      17_18     0.868    32.08 8.9e-05   5.67
3224           RPA2       1_19     0.867    21.72 6.0e-05   4.19
5262         SLAIN1      13_38     0.841    20.81 5.6e-05  -4.04
8242           NRG4      15_35     0.839    25.48 6.8e-05   4.48
10505       UGT2B17       4_48     0.835    48.94 1.3e-04   8.27
562           DGAT2      11_42     0.825   113.88 3.0e-04  13.70

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 29303.85 9.4e-02 -4.89
5520      RCN3      19_34         0  9358.83 0.0e+00 -5.07
8926   RPS6KB2      11_37         0  7563.91 2.9e-07 -4.55
7978    NDUFV1      11_37         0  6981.82 1.3e-12  1.53
168      SPRTN      1_118         0  5167.34 0.0e+00 -3.22
1818     ESRP2      16_36         0  4689.55 0.0e+00 -1.93
10901   PSMB10      16_36         0  4457.24 0.0e+00 -1.88
1794     NUTF2      16_36         0  4455.38 0.0e+00  2.04
6809  C16orf86      16_36         0  4399.66 0.0e+00 -1.89
1805       ACD      16_36         0  4363.60 0.0e+00 -1.92
1804      CTCF      16_36         0  3755.01 0.0e+00  1.73
374       EDC4      16_36         0  3728.73 0.0e+00  1.10
3138     EXOC8      1_118         0  3716.67 0.0e+00 -3.02
1796     CENPT      16_36         0  3691.79 0.0e+00  1.21
806     NFATC3      16_36         0  3642.73 0.0e+00 -2.76
6806  ATP6V0D1      16_36         0  3486.84 0.0e+00 -1.58
6805    ZDHHC1      16_36         0  3472.98 0.0e+00  1.48
7875      DUS2      16_36         0  3468.33 0.0e+00  6.75
6804     TPPP3      16_36         0  3376.47 0.0e+00  0.42
10904     E2F4      16_36         0  3318.89 0.0e+00  1.75

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 29303.85 0.09400  -4.89
5389           CTRL      16_36     0.979   258.47 0.00081  15.96
4564          PSRC1       1_67     0.998   231.34 0.00074  16.64
1267         PABPC4       1_24     0.868   185.89 0.00051  14.60
7089           USP1       1_39     0.550   255.43 0.00045  16.75
6089          FADS1      11_34     0.910   156.08 0.00045 -14.88
7654          PSMC3      11_29     0.997   109.99 0.00035 -14.44
4610           ACP2      11_29     0.996   105.08 0.00033 -14.33
3378           GPAM      10_70     0.973    99.00 0.00031  10.46
562           DGAT2      11_42     0.825   113.88 0.00030  13.70
7786       CATSPER2      15_16     0.669   116.07 0.00025 -10.07
10765       ZDHHC18       1_18     1.000    76.68 0.00024 -10.10
2442          IFT20      17_18     0.873    86.29 0.00024  -9.35
1652          PCIF1      20_28     0.999    64.50 0.00021   8.15
9322             F2      11_28     0.898    65.51 0.00019 -10.46
6959        CCDC116       22_4     0.492   106.26 0.00017  10.44
12237 RP11-110I1.14      11_71     0.994    54.36 0.00017   8.05
2496           ZPR1      11_70     0.140   364.62 0.00016   3.93
11441         APOC2      19_31     0.319   160.14 0.00016 -11.79
2410            MLX      17_25     0.882    55.30 0.00016  -7.32

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
8899       LPL       8_21     0.000  457.84 2.3e-13  30.54
3237     APOA1      11_70     0.000  197.52 0.0e+00 -17.89
7089      USP1       1_39     0.550  255.43 4.5e-04  16.75
4564     PSRC1       1_67     0.998  231.34 7.4e-04  16.64
3102     DOCK7       1_39     0.057  240.10 4.4e-05  16.05
5389      CTRL      16_36     0.979  258.47 8.1e-04  15.96
11020     LCAT      16_36     0.000  251.55 5.6e-09  15.35
6089     FADS1      11_34     0.910  156.08 4.5e-04 -14.88
1267    PABPC4       1_24     0.868  185.89 5.1e-04  14.60
5390     DPEP3      16_36     0.000  172.87 1.1e-15 -14.55
7654     PSMC3      11_29     0.997  109.99 3.5e-04 -14.44
6813     PSKH1      16_36     0.000 1096.57 1.0e-08 -14.36
4610      ACP2      11_29     0.996  105.08 3.3e-04 -14.33
4068   ALDH1A2      15_27     0.000  240.92 0.0e+00  13.99
6808   CARMIL2      16_36     0.000  200.07 1.5e-16 -13.82
562      DGAT2      11_42     0.825  113.88 3.0e-04  13.70
7874     DPEP2      16_36     0.000  198.04 2.9e-18  13.54
1807    PARD6A      16_36     0.000  171.75 3.8e-18 -13.47
7671      LIPC      15_27     0.000  535.17 1.3e-13  13.28
5500      AKT1      14_55     0.003  169.98 1.6e-06 -13.07

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.029653
#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
8899       LPL       8_21     0.000  457.84 2.3e-13  30.54
3237     APOA1      11_70     0.000  197.52 0.0e+00 -17.89
7089      USP1       1_39     0.550  255.43 4.5e-04  16.75
4564     PSRC1       1_67     0.998  231.34 7.4e-04  16.64
3102     DOCK7       1_39     0.057  240.10 4.4e-05  16.05
5389      CTRL      16_36     0.979  258.47 8.1e-04  15.96
11020     LCAT      16_36     0.000  251.55 5.6e-09  15.35
6089     FADS1      11_34     0.910  156.08 4.5e-04 -14.88
1267    PABPC4       1_24     0.868  185.89 5.1e-04  14.60
5390     DPEP3      16_36     0.000  172.87 1.1e-15 -14.55
7654     PSMC3      11_29     0.997  109.99 3.5e-04 -14.44
6813     PSKH1      16_36     0.000 1096.57 1.0e-08 -14.36
4610      ACP2      11_29     0.996  105.08 3.3e-04 -14.33
4068   ALDH1A2      15_27     0.000  240.92 0.0e+00  13.99
6808   CARMIL2      16_36     0.000  200.07 1.5e-16 -13.82
562      DGAT2      11_42     0.825  113.88 3.0e-04  13.70
7874     DPEP2      16_36     0.000  198.04 2.9e-18  13.54
1807    PARD6A      16_36     0.000  171.75 3.8e-18 -13.47
7671      LIPC      15_27     0.000  535.17 1.3e-13  13.28
5500      AKT1      14_55     0.003  169.98 1.6e-06 -13.07

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: 8_21"
       genename region_tag susie_pip    mu2     PVE     z
5936 CSGALNACT1       8_21         0  43.59 0.0e+00  6.66
1947     INTS10       8_21         0 133.94 0.0e+00  7.40
8899        LPL       8_21         0 457.84 2.3e-13 30.54

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 11_70"
     genename region_tag susie_pip    mu2     PVE      z
5007    BUD13      11_70      0.00 194.01 0.00000   1.66
2496     ZPR1      11_70      0.14 364.62 0.00016   3.93
3237    APOA1      11_70      0.00 197.52 0.00000 -17.89
6898     SIK3      11_70      0.00  38.39 0.00000  -5.62
8030 PAFAH1B2      11_70      0.00  82.61 0.00000  -4.60
6104    TAGLN      11_70      0.00 356.64 0.00000  -4.49
6902    PCSK7      11_70      0.00 287.15 0.00000  11.90
7873   RNF214      11_70      0.00  10.17 0.00000  -0.23
9915    BACE1      11_70      0.00  45.14 0.00000   0.28
2530   CEP164      11_70      0.00  74.71 0.00000   4.21
5018    FXYD2      11_70      0.00  17.75 0.00000   1.77
5017    FXYD6      11_70      0.00  12.87 0.00000   1.86

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 1_39"
     genename region_tag susie_pip    mu2     PVE     z
7088    TM2D1       1_39     0.023   4.95 3.7e-07  0.35
4449     PATJ       1_39     0.097  19.23 6.0e-06 -2.17
7089     USP1       1_39     0.550 255.43 4.5e-04 16.75
3102    DOCK7       1_39     0.057 240.10 4.4e-05 16.05
3822    ATG4C       1_39     0.067  16.63 3.6e-06 -2.34

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 1_67"
          genename region_tag susie_pip    mu2     PVE     z
11280 RP11-356N1.2       1_67     0.011   6.10 2.2e-07  0.95
1102      SLC25A24       1_67     0.022  11.36 8.1e-07 -1.24
7095       FAM102B       1_67     0.013   6.18 2.5e-07  0.54
7096        HENMT1       1_67     0.131  26.71 1.1e-05 -2.90
3080        STXBP3       1_67     0.027  15.56 1.3e-06 -2.08
3522         GPSM2       1_67     0.020  10.51 6.8e-07  1.24
3521         CLCC1       1_67     0.012   6.06 2.2e-07  0.15
10487        TAF13       1_67     0.015   8.85 4.2e-07  1.86
11143     TMEM167B       1_67     0.153  33.42 1.6e-05 -3.93
9291      C1orf194       1_67     0.019  12.03 7.1e-07  1.85
1099         WDR47       1_67     0.017  11.21 6.0e-07  1.80
3084      KIAA1324       1_67     0.012  10.17 3.9e-07 -2.10
331           SARS       1_67     0.020  33.72 2.2e-06  5.10
5562        CELSR2       1_67     0.044  24.38 3.4e-06 -4.39
4564         PSRC1       1_67     0.998 231.34 7.4e-04 16.64
7099       ATXN7L2       1_67     0.013   6.56 2.7e-07 -0.30
8776      CYB561D1       1_67     0.018  13.27 7.4e-07  1.56
9435        AMIGO1       1_67     0.028  13.20 1.2e-06  0.74
617          GNAI3       1_67     0.091  31.99 9.3e-06 -3.42
11016        GSTM2       1_67     0.017  50.60 2.7e-06 -7.02
8107         GSTM4       1_67     0.020  11.27 7.1e-07  0.55
4559         GSTM1       1_67     0.011  61.75 2.1e-06 -8.92
4561         GSTM5       1_67     0.011  13.67 4.7e-07 -3.83
4562         GSTM3       1_67     0.016   9.88 5.0e-07 -0.12

Version Author Date
dfd2b5f wesleycrouse 2021-09-07
[1] "Region: 16_36"
           genename region_tag susie_pip     mu2     PVE      z
4752       DYNC1LI2      16_36     0.000   23.54 0.0e+00   1.06
9282           CDH5      16_36     0.000   31.02 0.0e+00  -2.97
11678     LINC00920      16_36     0.000   78.37 0.0e+00  -1.59
7763          BEAN1      16_36     0.000   35.42 0.0e+00   3.40
7764            TK2      16_36     0.000   24.58 0.0e+00   2.86
11156          CKLF      16_36     0.000   59.70 0.0e+00   2.43
1233          CMTM1      16_36     0.000   36.19 0.0e+00   3.39
5365          CMTM3      16_36     0.000  131.44 0.0e+00   3.78
9636          CMTM4      16_36     0.000   22.78 0.0e+00   3.21
6794           NAE1      16_36     0.000   21.70 0.0e+00  -1.71
8627           PDP2      16_36     0.000   28.15 0.0e+00  -0.52
8626           CES2      16_36     0.000  163.62 0.0e+00  -2.36
8624           CES3      16_36     0.000   24.15 0.0e+00   2.66
695            CBFB      16_36     0.000   44.55 0.0e+00  -5.52
3773       C16orf70      16_36     0.000   43.14 0.0e+00  -5.47
11479        B3GNT9      16_36     0.000  157.31 0.0e+00   1.22
5366           NOL3      16_36     0.000   45.73 0.0e+00   4.63
1793          ELMO3      16_36     0.000   64.09 0.0e+00   6.78
10210     KIAA0895L      16_36     0.000    8.74 0.0e+00   0.08
9235        EXOC3L1      16_36     0.000   60.61 0.0e+00   6.69
10904          E2F4      16_36     0.000 3318.89 0.0e+00   1.75
4756         SLC9A5      16_36     0.000   51.32 0.0e+00   5.90
3769         LRRC29      16_36     0.000   66.26 0.0e+00  -6.87
4754          FHOD1      16_36     0.000   17.88 0.0e+00   0.89
10218       PLEKHG4      16_36     0.000 1832.66 0.0e+00   4.88
6804          TPPP3      16_36     0.000 3376.47 0.0e+00   0.42
6805         ZDHHC1      16_36     0.000 3472.98 0.0e+00   1.48
6806       ATP6V0D1      16_36     0.000 3486.84 0.0e+00  -1.58
1804           CTCF      16_36     0.000 3755.01 0.0e+00   1.73
12029 CTD-2012K14.6      16_36     0.000   49.51 0.0e+00  -0.07
6808        CARMIL2      16_36     0.000  200.07 1.5e-16 -13.82
1807         PARD6A      16_36     0.000  171.75 3.8e-18 -13.47
1805            ACD      16_36     0.000 4363.60 0.0e+00  -1.92
3665          ENKD1      16_36     0.000   30.49 0.0e+00  -1.05
6809       C16orf86      16_36     0.000 4399.66 0.0e+00  -1.89
5391          GFOD2      16_36     0.000   18.72 0.0e+00   0.72
1797       TSNAXIP1      16_36     0.000   19.42 0.0e+00   1.67
1794          NUTF2      16_36     0.000 4455.38 0.0e+00   2.04
1796          CENPT      16_36     0.000 3691.79 0.0e+00   1.21
374            EDC4      16_36     0.000 3728.73 0.0e+00   1.10
10064         NRN1L      16_36     0.000   73.84 0.0e+00  -3.85
6813          PSKH1      16_36     0.000 1096.57 1.0e-08 -14.36
5389           CTRL      16_36     0.979  258.47 8.1e-04  15.96
10901        PSMB10      16_36     0.000 4457.24 0.0e+00  -1.88
5390          DPEP3      16_36     0.000  172.87 1.1e-15 -14.55
11020          LCAT      16_36     0.000  251.55 5.6e-09  15.35
7875           DUS2      16_36     0.000 3468.33 0.0e+00   6.75
7874          DPEP2      16_36     0.000  198.04 2.9e-18  13.54
806          NFATC3      16_36     0.000 3642.73 0.0e+00  -2.76
1818          ESRP2      16_36     0.000 4689.55 0.0e+00  -1.93
1816         SLC7A6      16_36     0.000  174.59 0.0e+00 -12.59
1817        PLA2G15      16_36     0.000  128.26 0.0e+00 -10.66
4414          PRMT7      16_36     0.000   75.15 0.0e+00   8.98
9744          ZFP90      16_36     0.000   21.50 0.0e+00   1.04

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
32910     rs9427104       1_75     1.000    59.94 1.9e-04   7.67
52630     rs2642420      1_112     1.000    42.91 1.4e-04  -7.64
55284     rs2103827      1_117     1.000   177.59 5.7e-04  19.40
55285    rs11122453      1_117     1.000   355.61 1.1e-03  22.94
55766   rs766167074      1_118     1.000  5284.70 1.7e-02   3.00
67964     rs1042034       2_13     1.000   439.62 1.4e-03 -22.04
67965     rs1801699       2_13     1.000    40.04 1.3e-04  -7.66
69700      rs780093       2_16     1.000   107.83 3.4e-04 -11.04
179579    rs9817452       3_97     1.000    53.13 1.7e-04   7.45
226226   rs35518360       4_67     1.000   218.42 7.0e-04 -15.55
226292   rs13140033       4_68     1.000   127.17 4.1e-04 -11.78
270843   rs62369502       5_28     1.000    47.70 1.5e-04  -6.90
365965  rs191555775      6_104     1.000    61.69 2.0e-04  -9.19
411549    rs6961342       7_80     1.000   149.04 4.8e-04 -11.19
426292    rs7012814       8_12     1.000   126.40 4.0e-04  13.90
426595    rs1402522       8_13     1.000    52.90 1.7e-04   7.83
427356   rs17810889       8_15     1.000    54.96 1.8e-04   7.91
431557   rs75835816       8_21     1.000   432.77 1.4e-03 -21.65
431593   rs11986461       8_21     1.000   297.44 9.5e-04  19.74
462183   rs10956254       8_83     1.000    34.53 1.1e-04  -6.15
468461    rs7832515       8_94     1.000   145.72 4.6e-04  12.84
495023    rs2777798       9_52     1.000   372.16 1.2e-03  14.36
495031    rs2777804       9_52     1.000   144.15 4.6e-04   4.64
495037    rs7024300       9_53     1.000   253.55 8.1e-04  16.99
495043    rs2297400       9_53     1.000   256.21 8.2e-04  15.59
495054   rs62568181       9_53     1.000   339.56 1.1e-03 -24.86
495067    rs2254819       9_53     1.000   252.37 8.1e-04 -22.66
495070    rs2437818       9_53     1.000   112.25 3.6e-04  15.45
518864   rs71007692      10_28     1.000  4974.96 1.6e-02  -2.78
578713     rs640621      11_70     1.000   781.61 2.5e-03 -24.74
601920    rs6581124      12_35     1.000    52.68 1.7e-04   8.14
601939    rs7397189      12_36     1.000    66.46 2.1e-04  11.12
605980    rs2137537      12_44     1.000    40.50 1.3e-04  -5.98
622624   rs61941677      12_76     1.000   140.43 4.5e-04 -14.41
622630   rs12230146      12_76     1.000    60.74 1.9e-04   7.50
674672   rs13379043      14_34     1.000    67.67 2.2e-04   8.46
695142   rs72737411      15_25     1.000    43.28 1.4e-04  -6.21
695252   rs58038553      15_27     1.000   342.57 1.1e-03 -24.93
695254    rs1711037      15_27     1.000   153.87 4.9e-04  16.41
695312   rs28594460      15_27     1.000   334.39 1.1e-03  21.14
695328   rs62000868      15_27     1.000   894.23 2.9e-03  31.87
695334    rs2070895      15_27     1.000  2223.91 7.1e-03  50.25
721595    rs8064102      16_31     1.000   294.95 9.4e-04   6.16
721617  rs190575415      16_31     1.000   311.73 9.9e-04  16.53
721627     rs821840      16_31     1.000  4452.45 1.4e-02  74.05
721628   rs12720926      16_31     1.000  2773.56 8.9e-03  66.58
721632   rs66495554      16_31     1.000   985.37 3.1e-03  -7.44
724451    rs2276329      16_37     1.000    49.88 1.6e-04  -6.87
754781  rs146424771       18_3     1.000    64.89 2.1e-04   2.47
767463   rs11082766      18_27     1.000   181.93 5.8e-04  13.66
767483    rs6507938      18_27     1.000   666.13 2.1e-03  32.32
767484  rs118043171      18_27     1.000   412.24 1.3e-03  27.52
767703   rs74461650      18_28     1.000   109.21 3.5e-04  10.92
782307    rs1865063      19_11     1.000   210.52 6.7e-04 -14.13
782309    rs3745683      19_11     1.000    89.99 2.9e-04 -13.33
791435     rs814573      19_31     1.000   689.64 2.2e-03 -29.76
791441    rs4803775      19_31     1.000   186.72 6.0e-04  14.22
806860    rs6063139      20_29     1.000    75.34 2.4e-04   3.79
806898   rs78492788      20_29     1.000    94.99 3.0e-04   9.49
872966  rs140584594       1_67     1.000   159.28 5.1e-04  14.42
920095   rs35733538       3_95     1.000  2313.38 7.4e-03  -5.36
926460   rs35945848        4_2     1.000  1953.28 6.2e-03   3.13
986221   rs11601507       11_4     1.000    51.94 1.7e-04  -6.79
1023190 rs146923372      11_37     1.000 14480.70 4.6e-02   4.02
1084266   rs4986970      16_36     1.000   145.96 4.7e-04 -12.74
1085097  rs56090907      16_36     1.000  6377.25 2.0e-02   4.77
1116082  rs72836561      17_26     1.000   444.60 1.4e-03 -21.66
1119239 rs116843064       19_8     1.000   225.35 7.2e-04  16.71
1137782  rs61371437      19_34     1.000 28186.41 9.0e-02   4.89
1137794 rs374141296      19_34     1.000 28260.93 9.0e-02   5.19
1152147  rs12975366      19_37     1.000    79.79 2.5e-04 -10.99
1163659   rs1800961      20_28     1.000   441.16 1.4e-03 -21.94
325282  rs181268076       6_27     0.999    61.80 2.0e-04  -7.07
417541    rs6977416       7_94     0.999    39.35 1.3e-04  -6.32
426303   rs13265179       8_12     0.999   162.85 5.2e-04 -16.86
622608    rs3782287      12_76     0.999    62.45 2.0e-04 -11.72
622638   rs11834751      12_76     0.999    48.86 1.6e-04  -4.27
728789   rs12443634      16_46     0.999    69.01 2.2e-04  10.14
791402  rs111794050      19_31     0.999    70.30 2.2e-04   9.35
791431     rs405509      19_31     0.999    78.58 2.5e-04  13.56
791759   rs58701309      19_32     0.999    86.20 2.7e-04   3.46
1023185  rs57808037      11_37     0.999 14479.15 4.6e-02   4.04
782338   rs35753011      19_11     0.998    84.91 2.7e-04  -0.92
5647       rs603412       1_14     0.997    38.30 1.2e-04  -6.26
31792   rs185073199       1_73     0.997    31.10 9.9e-05   5.49
426291     rs713286       8_12     0.997    54.78 1.7e-04 -12.24
503785  rs115478735       9_70     0.997    91.91 2.9e-04   9.87
460707    rs9297630       8_80     0.996    55.21 1.8e-04  -7.19
498390    rs2763193       9_59     0.996    51.42 1.6e-04   5.91
767499   rs62101781      18_27     0.996   299.90 9.5e-04  19.53
207416   rs58932203       4_32     0.995    43.64 1.4e-04   6.15
323700   rs77424484       6_26     0.994    66.81 2.1e-04  -7.74
748678  rs113408695      17_39     0.993    38.01 1.2e-04  -5.86
34934     rs4657041       1_79     0.992    32.76 1.0e-04  -5.77
325011    rs2858317       6_26     0.992   116.83 3.7e-04  11.56
543302   rs10901802      10_78     0.992    33.94 1.1e-04   5.87
653598  rs185932947      13_52     0.992    27.71 8.8e-05   5.10
790334   rs11879413      19_30     0.992    31.74 1.0e-04   5.64
61939    rs10183939        2_2     0.991    32.06 1.0e-04  -5.59
228047  rs111349657       4_71     0.991    30.53 9.7e-05  -5.71
498382    rs7040440       9_59     0.991    30.71 9.7e-05   3.18
1062725   rs2494732      14_55     0.991   102.70 3.2e-04   1.84
1115618 rs117380643      17_25     0.991    68.72 2.2e-04  -8.28
697062   rs11071771      15_29     0.989    35.61 1.1e-04  -5.86
328541  rs142449754       6_32     0.988    36.29 1.1e-04  -6.44
581366    rs4937122      11_77     0.988    32.22 1.0e-04  -5.50
625794   rs76734539      12_82     0.988    30.81 9.7e-05   5.68
579848    rs1219430      11_74     0.982    30.60 9.6e-05  -5.72
708180    rs1037118      15_50     0.982    28.59 9.0e-05   4.67
142353    rs4681065       3_19     0.981    28.79 9.0e-05   5.17
753084   rs72854483      17_46     0.976    26.25 8.2e-05  -4.92
701050   rs16972386      15_38     0.974    27.38 8.5e-05  -4.97
791760    rs7259871      19_32     0.974   100.61 3.1e-04  -5.78
169035     rs334563       3_74     0.973    46.18 1.4e-04   6.67
864778   rs11591147       1_34     0.972    31.71 9.8e-05   5.65
322856    rs3095311       6_26     0.970   115.58 3.6e-04 -10.67
370246   rs11971790        7_3     0.970    59.52 1.8e-04  -6.93
952940     rs686030       9_13     0.970   241.76 7.5e-04  16.56
384656       rs9490       7_28     0.969    26.92 8.3e-05   4.97
574302   rs72980276      11_59     0.967    26.70 8.2e-05  -4.99
556662   rs12288512      11_19     0.964    53.72 1.7e-04  -7.27
323880    rs3869145       6_26     0.963    55.39 1.7e-04  -8.02
67911    rs17721572       2_12     0.962    36.45 1.1e-04  -4.63
380411   rs75348547       7_22     0.962    27.65 8.5e-05  -4.86
478373  rs145804707       9_18     0.960    24.96 7.6e-05  -4.69
728796   rs11641142      16_46     0.960    48.51 1.5e-04   9.01
321497    rs1233480       6_23     0.958    94.78 2.9e-04  -9.65
179337     rs816547       3_97     0.953    25.37 7.7e-05  -4.76
272758  rs113088001       5_31     0.953    32.70 9.9e-05  -5.21
431539    rs6999813       8_21     0.950   435.42 1.3e-03  31.97
812524  rs759884584      20_38     0.942    31.03 9.3e-05   2.74
52628      rs883847      1_112     0.941    33.59 1.0e-04   6.84
193684   rs13116176        4_4     0.940    29.38 8.8e-05  -4.46
601875    rs1874888      12_35     0.940    33.30 1.0e-04  -6.19
380588    rs2699814       7_23     0.934    27.09 8.1e-05   4.98
273694     rs173964       5_33     0.928    80.30 2.4e-04  -7.90
995001    rs7123635      11_28     0.927    47.77 1.4e-04  -7.04
69407   rs577641882       2_15     0.926    33.03 9.8e-05  -5.71
207327   rs12639940       4_32     0.926    24.19 7.2e-05  -3.63
770359   rs41292412      18_31     0.926    30.52 9.0e-05  -5.49
301991    rs4958365       5_90     0.925    31.09 9.2e-05   4.70
321253    rs3132390       6_22     0.919    86.91 2.5e-04  -9.77
431331  rs113231830       8_20     0.919    26.10 7.6e-05  -4.89
53794    rs12132342      1_115     0.918    24.45 7.2e-05  -4.53
132476    rs4675812      2_144     0.910    32.67 9.5e-05   6.69
329293     rs581465       6_34     0.909    26.03 7.5e-05   4.83
1021356   rs4930352      11_37     0.907   461.61 1.3e-03   7.80
741680   rs12450388      17_23     0.904    25.34 7.3e-05   4.39
362643    rs6905582       6_99     0.902    24.27 7.0e-05  -4.49
828336      rs12321       22_9     0.900    25.09 7.2e-05   4.48
525850   rs10761739      10_42     0.895    41.45 1.2e-04   6.40
590855     rs964974      12_15     0.892    36.91 1.1e-04  -5.99
622601     rs838876      12_76     0.889   126.44 3.6e-04 -13.65
832151  rs531420711      22_17     0.886    25.45 7.2e-05  -4.81
606101    rs1707498      12_44     0.882    29.31 8.2e-05   5.11
355180  rs112120898       6_84     0.869    24.37 6.8e-05   4.42
431523   rs17091881       8_21     0.866   315.07 8.7e-04 -18.25
550501    rs2923096       11_8     0.864    33.26 9.2e-05   5.90
459205    rs2721932       8_78     0.863   248.86 6.9e-04  14.87
578686    rs3135506      11_70     0.860   357.80 9.8e-04  -4.60
586800   rs10849492       12_7     0.856    29.13 8.0e-05  -5.34
793518   rs34637812      19_38     0.856    25.67 7.0e-05  -4.61
284264  rs538659275       5_57     0.855    26.01 7.1e-05   4.87
459259   rs10095930       8_78     0.855    75.60 2.1e-04   5.62
601962  rs140734681      12_36     0.855    22.23 6.1e-05  -1.53
55277     rs6678475      1_117     0.853    25.93 7.1e-05  -1.21
908148  rs115787626       3_33     0.853    54.64 1.5e-04  -7.73
439582   rs62515418       8_40     0.850    24.57 6.7e-05   4.31
673133     rs194749      14_33     0.849    41.11 1.1e-04   6.09
155733   rs73090625       3_48     0.847    30.56 8.3e-05   5.31
685827   rs11161243       15_4     0.847    32.36 8.7e-05   5.48
294648   rs10057561       5_77     0.846    26.73 7.2e-05  -4.91
1169161    rs235314      21_23     0.846    85.37 2.3e-04  -9.34
84728     rs3796098       2_47     0.843    30.28 8.1e-05  -5.43
149031   rs78342753       3_35     0.842    33.26 8.9e-05  -6.36
794824    rs2316866       20_1     0.839    25.29 6.8e-05  -4.54
325153    rs9276189       6_27     0.835    37.27 9.9e-05   3.68
541420   rs78787582      10_74     0.833    25.93 6.9e-05   4.69
517323  rs145407036      10_24     0.829    29.01 7.7e-05  -5.44
738786   rs78792395      17_15     0.828    26.33 7.0e-05   4.54
71055    rs11124265       2_20     0.817    25.39 6.6e-05   4.48
455512    rs2570952       8_69     0.813    30.32 7.9e-05   5.25
468451    rs7387969       8_94     0.813    40.17 1.0e-04   6.84
368701   rs34207042      6_110     0.812    24.70 6.4e-05   4.28
211185     rs723585       4_40     0.811    36.57 9.5e-05   5.82
767329   rs76055069      18_27     0.811   319.03 8.3e-04  25.71
517546   rs12765967      10_25     0.807    29.61 7.6e-05  -5.17
202380   rs56147366       4_22     0.804    44.29 1.1e-04  -6.56
460843    rs4871180       8_80     0.802    26.93 6.9e-05   4.24

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
1137794 rs374141296      19_34         1 28260.93 9.0e-02 5.19
1137791 rs113176985      19_34         0 28214.13 0.0e+00 4.85
1137798   rs2946865      19_34         0 28202.09 0.0e+00 4.89
1137782  rs61371437      19_34         1 28186.41 9.0e-02 4.89
1137784  rs35295508      19_34         0 28181.31 0.0e+00 4.89
1137789  rs73056069      19_34         0 28079.60 0.0e+00 4.87
1137773    rs756628      19_34         0 28071.84 2.8e-07 4.88
1137772    rs739349      19_34         0 28071.67 2.0e-07 4.88
1137769    rs739347      19_34         0 28011.77 9.6e-14 4.79
1137786   rs2878354      19_34         0 28009.82 0.0e+00 4.91
1137770   rs2073614      19_34         0 27974.95 3.0e-17 4.75
1137775   rs2077300      19_34         0 27903.13 1.2e-13 4.90
1137765   rs4802613      19_34         0 27851.90 0.0e+00 4.75
1137779  rs73056059      19_34         0 27850.11 0.0e+00 4.84
1137799  rs60815603      19_34         0 27709.16 0.0e+00 5.04
1137802   rs1316885      19_34         0 27673.86 0.0e+00 5.13
1137807   rs2946863      19_34         0 27625.22 0.0e+00 5.19
1137800  rs35443645      19_34         0 27580.12 0.0e+00 5.05
1137804  rs60746284      19_34         0 27538.77 0.0e+00 5.17
1137763  rs10403394      19_34         0 27482.16 0.0e+00 4.94
1137764  rs17555056      19_34         0 27459.82 0.0e+00 4.87
1137780  rs73056062      19_34         0 27131.74 0.0e+00 4.80
1137810 rs553431297      19_34         0 26738.29 0.0e+00 4.57
1137793 rs112283514      19_34         0 26649.08 0.0e+00 4.82
1137795  rs11270139      19_34         0 26481.17 0.0e+00 4.44
1137760  rs10421294      19_34         0 24805.48 0.0e+00 3.94
1137759   rs8108175      19_34         0 24801.93 0.0e+00 3.93
1137752  rs59192944      19_34         0 24754.20 0.0e+00 3.91
1137758   rs1858742      19_34         0 24743.70 0.0e+00 3.86
1137749  rs55991145      19_34         0 24733.42 0.0e+00 3.92
1137744   rs3786567      19_34         0 24723.57 0.0e+00 3.92
1137740   rs2271952      19_34         0 24713.32 0.0e+00 3.91
1137743   rs4801801      19_34         0 24709.97 0.0e+00 3.90
1137739   rs2271953      19_34         0 24681.35 0.0e+00 3.87
1137741   rs2271951      19_34         0 24680.11 0.0e+00 3.87
1137730  rs60365978      19_34         0 24662.80 0.0e+00 3.90
1137736   rs4802612      19_34         0 24563.69 0.0e+00 3.99
1137746   rs2517977      19_34         0 24494.38 0.0e+00 3.81
1137733  rs55893003      19_34         0 24483.30 0.0e+00 3.93
1137725  rs55992104      19_34         0 23917.19 0.0e+00 3.75
1137719  rs60403475      19_34         0 23915.37 0.0e+00 3.77
1137722   rs4352151      19_34         0 23906.18 0.0e+00 3.73
1137716  rs11878448      19_34         0 23888.90 0.0e+00 3.72
1137710   rs9653100      19_34         0 23884.66 0.0e+00 3.74
1137706   rs4802611      19_34         0 23870.19 0.0e+00 3.77
1137698   rs7251338      19_34         0 23833.66 0.0e+00 3.76
1137697  rs59269605      19_34         0 23831.12 0.0e+00 3.77
1137718   rs1042120      19_34         0 23766.04 0.0e+00 3.82
1137714 rs113220577      19_34         0 23744.64 0.0e+00 3.81
1137708   rs9653118      19_34         0 23706.52 0.0e+00 3.78

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
1137782  rs61371437      19_34     1.000 28186.41 0.0900   4.89
1137794 rs374141296      19_34     1.000 28260.93 0.0900   5.19
1023185  rs57808037      11_37     0.999 14479.15 0.0460   4.04
1023190 rs146923372      11_37     1.000 14480.70 0.0460   4.02
1085097  rs56090907      16_36     1.000  6377.25 0.0200   4.77
55766   rs766167074      1_118     1.000  5284.70 0.0170   3.00
518864   rs71007692      10_28     1.000  4974.96 0.0160  -2.78
721627     rs821840      16_31     1.000  4452.45 0.0140  74.05
721628   rs12720926      16_31     1.000  2773.56 0.0089  66.58
1085067  rs71395853      16_36     0.436  6402.29 0.0089   1.49
518870    rs2472183      10_28     0.482  5005.57 0.0077  -2.84
518863    rs2474565      10_28     0.469  5005.53 0.0075  -2.84
920095   rs35733538       3_95     1.000  2313.38 0.0074  -5.36
518873   rs11011452      10_28     0.458  5005.63 0.0073  -2.83
695334    rs2070895      15_27     1.000  2223.91 0.0071  50.25
926460   rs35945848        4_2     1.000  1953.28 0.0062   3.13
518861    rs9299760      10_28     0.335  5003.06 0.0053  -2.85
55763    rs10489611      1_118     0.298  5324.45 0.0051   3.25
1085099  rs71395854      16_36     0.237  6403.15 0.0048   1.46
55765      rs971534      1_118     0.270  5324.37 0.0046   3.25
55757     rs2256908      1_118     0.256  5324.16 0.0043   3.25
55764     rs2486737      1_118     0.248  5324.32 0.0042   3.24
55753     rs1076804      1_118     0.218  5317.93 0.0037   3.29
55760     rs2790891      1_118     0.196  5324.00 0.0033   3.24
55761     rs2491405      1_118     0.196  5324.00 0.0033   3.24
721632   rs66495554      16_31     1.000   985.37 0.0031  -7.44
1085110  rs35189054      16_36     0.154  6402.06 0.0031   1.46
695328   rs62000868      15_27     1.000   894.23 0.0029  31.87
578713     rs640621      11_70     1.000   781.61 0.0025 -24.74
926479    rs1680073        4_2     0.383  1948.46 0.0024   2.76
791435     rs814573      19_31     1.000   689.64 0.0022 -29.76
767483    rs6507938      18_27     1.000   666.13 0.0021  32.32
926470    rs1680074        4_2     0.344  1948.55 0.0021   2.75
1023200   rs6591245      11_37     0.043 14456.05 0.0020   3.96
67964     rs1042034       2_13     1.000   439.62 0.0014 -22.04
431557   rs75835816       8_21     1.000   432.77 0.0014 -21.65
920099     rs433903       3_95     0.196  2242.10 0.0014  -5.44
920104     rs355782       3_95     0.194  2242.34 0.0014  -5.44
926486   rs10024013        4_2     0.218  1948.01 0.0014   2.72
1085069  rs34530665      16_36     0.067  6401.82 0.0014   1.45
1116082  rs72836561      17_26     1.000   444.60 0.0014 -21.66
1163659   rs1800961      20_28     1.000   441.16 0.0014 -21.94
431539    rs6999813       8_21     0.950   435.42 0.0013  31.97
767484  rs118043171      18_27     1.000   412.24 0.0013  27.52
1021356   rs4930352      11_37     0.907   461.61 0.0013   7.80
495023    rs2777798       9_52     1.000   372.16 0.0012  14.36
55285    rs11122453      1_117     1.000   355.61 0.0011  22.94
495054   rs62568181       9_53     1.000   339.56 0.0011 -24.86
695252   rs58038553      15_27     1.000   342.57 0.0011 -24.93
695312   rs28594460      15_27     1.000   334.39 0.0011  21.14

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
721627    rs821840      16_31     1.000 4452.45 1.4e-02  74.05
721625  rs12446515      16_31     0.000 4410.80 5.9e-12  73.66
721628  rs12720926      16_31     1.000 2773.56 8.9e-03  66.58
721624    rs193695      16_31     0.000 2028.32 8.7e-15  50.33
695334   rs2070895      15_27     1.000 2223.91 7.1e-03  50.25
695336   rs8034802      15_27     0.000 1439.98 0.0e+00  39.84
695347    rs686958      15_27     0.000 1387.04 0.0e+00 -39.62
695337    rs633695      15_27     0.000 1365.45 0.0e+00  38.77
695340    rs261341      15_27     0.000 1309.73 0.0e+00 -37.82
695342    rs488490      15_27     0.000 1335.99 0.0e+00 -37.73
695346    rs485671      15_27     0.000 1305.03 0.0e+00 -37.05
431534   rs2410620       8_21     0.469  428.07 6.4e-04  35.48
431541   rs1441762       8_21     0.311  426.37 4.2e-04  35.46
431544   rs4126104       8_21     0.040  420.58 5.3e-05  35.41
431540  rs35878331       8_21     0.124  422.11 1.7e-04  35.35
431545  rs35369244       8_21     0.004  412.71 5.2e-06  35.09
431514 rs149553676       8_21     0.000  311.99 0.0e+00  32.53
431513       rs287       8_21     0.000  317.58 0.0e+00  32.42
767483   rs6507938      18_27     1.000  666.13 2.1e-03  32.32
431543  rs11986942       8_21     0.052  457.17 7.6e-05  32.15
431518  rs10645926       8_21     0.000  435.68 6.0e-07  31.99
431539   rs6999813       8_21     0.950  435.42 1.3e-03  31.97
695328  rs62000868      15_27     1.000  894.23 2.9e-03  31.87
767480   rs7241918      18_27     0.000  640.63 1.6e-08  31.84
431525  rs78963197       8_21     0.048  434.09 6.7e-05  31.58
431530  rs17489226       8_21     0.001  415.56 8.1e-07  31.38
431551  rs28675909       8_21     0.000  412.58 2.2e-07  31.31
431554  rs11989309       8_21     0.000  411.78 1.6e-07  31.30
431556   rs7816447       8_21     0.000  415.44 3.0e-07  31.30
431552  rs79198716       8_21     0.000  409.14 8.8e-08  31.27
431553   rs7004149       8_21     0.000  408.98 8.3e-08  31.26
431558  rs11984698       8_21     0.000  409.18 7.1e-08  31.24
578802   rs7930783      11_70     0.330  631.10 6.7e-04  31.22
578804   rs7932655      11_70     0.328  631.03 6.6e-04  31.22
578807  rs59097294      11_70     0.332  631.04 6.7e-04  31.22
431520   rs1569209       8_21     0.000  397.73 8.2e-14  30.70
431522  rs80073370       8_21     0.000  389.17 1.7e-15  30.35
578715  rs11216162      11_70     0.010  675.74 2.2e-05  30.18
431583  rs80026582       8_21     0.000  418.63 6.5e-16  30.05
791435    rs814573      19_31     1.000  689.64 2.2e-03 -29.76
721630   rs4784744      16_31     0.000 1194.42 0.0e+00 -29.10
721629    rs289717      16_31     0.000 1188.16 0.0e+00 -29.05
721631   rs4784745      16_31     0.000 1240.66 0.0e+00 -28.97
431536   rs4083261       8_21     0.000  406.06 1.3e-11  28.57
431535  rs12541912       8_21     0.000  362.71 8.9e-15  27.95
767466   rs4121823      18_27     0.000  548.60 0.0e+00  27.72
767484 rs118043171      18_27     1.000  412.24 1.3e-03  27.52
767329  rs76055069      18_27     0.811  319.03 8.3e-04  25.71
431548   rs6586886       8_21     0.000  166.81 0.0e+00  25.51
767326  rs75495141      18_27     0.189  310.96 1.9e-04  25.44

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] 40
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)
BLMH gene(s) from the input list not found in DisGeNET CURATEDMRPL21 gene(s) from the input list not found in DisGeNET CURATEDSLAIN1 gene(s) from the input list not found in DisGeNET CURATEDSP3 gene(s) from the input list not found in DisGeNET CURATEDPCIF1 gene(s) from the input list not found in DisGeNET CURATEDCEBPG gene(s) from the input list not found in DisGeNET CURATEDMARCH2 gene(s) from the input list not found in DisGeNET CURATEDFOXK1 gene(s) from the input list not found in DisGeNET CURATEDRPA2 gene(s) from the input list not found in DisGeNET CURATEDIFT20 gene(s) from the input list not found in DisGeNET CURATEDRP11-54O7.17 gene(s) from the input list not found in DisGeNET CURATEDNRG4 gene(s) from the input list not found in DisGeNET CURATEDNTAN1 gene(s) from the input list not found in DisGeNET CURATEDSLFN13 gene(s) from the input list not found in DisGeNET CURATEDFCGRT gene(s) from the input list not found in DisGeNET CURATEDLAMP1 gene(s) from the input list not found in DisGeNET CURATEDZDHHC18 gene(s) from the input list not found in DisGeNET CURATEDRP11-110I1.14 gene(s) from the input list not found in DisGeNET CURATEDSPSB1 gene(s) from the input list not found in DisGeNET CURATEDPSRC1 gene(s) from the input list not found in DisGeNET CURATED
                                    Description       FDR Ratio  BgRatio
27             Glycogen Storage Disease Type IV 0.0147245  1/20   1/9703
32                               HIV Infections 0.0147245  3/20 103/9703
34               Inherited Factor II deficiency 0.0147245  1/20   1/9703
67                      Skin Diseases, Vascular 0.0147245  1/20   1/9703
79                  Acid Phosphatase Deficiency 0.0147245  1/20   1/9703
80      Hereditary factor II deficiency disease 0.0147245  1/20   1/9703
118       Polyglucosan Body Disease, Adult Form 0.0147245  1/20   1/9703
119 GSD IV, Neuromuscular Form, Fatal Perinatal 0.0147245  1/20   1/9703
120      GSD IV, Neuromuscular Form, Congenital 0.0147245  1/20   1/9703
121       GSD IV, Neuromuscular Form, Childhood 0.0147245  1/20   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