Last updated: 2021-09-09
Checks: 7 0
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.
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
The command set.seed(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
Modified: analysis/ukb-d-30890_irnt_Whole_Blood.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-30660_irnt_Liver_trim.Rmd
) and HTML (docs/ukb-d-30660_irnt_Liver_trim.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 |
---|---|---|---|---|
html | cbf7408 | wesleycrouse | 2021-09-08 | adding enrichment to reports |
html | 4970e3e | wesleycrouse | 2021-09-08 | updating reports |
html | dfd2b5f | wesleycrouse | 2021-09-07 | regenerating reports |
html | 47f58ac | wesleycrouse | 2021-09-06 | fixing thin argument for fixed pi results |
Rmd | 209346f | wesleycrouse | 2021-09-06 | updating additional analyses |
html | 209346f | wesleycrouse | 2021-09-06 | updating additional analyses |
html | b14741c | wesleycrouse | 2021-09-06 | switching from render to wflow_build |
html | 61b53b3 | wesleycrouse | 2021-09-06 | updated PVE calculation |
html | 837dd01 | wesleycrouse | 2021-09-01 | adding additional fixedsigma report |
html | 2739f4f | wesleycrouse | 2021-08-30 | fixing typo |
html | b1e6b7e | wesleycrouse | 2021-08-30 | fixing alignment on index |
html | d7dfe76 | wesleycrouse | 2021-08-30 | Adding detailed reports for 30660 |
Rmd | ea2e654 | wesleycrouse | 2021-08-30 | Exploring fixed priors and trimming large z scores |
These are the results of a ctwas
analysis of the UK Biobank trait Direct bilirubin (quantile)
using Liver
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-30660_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 Liver
eQTL obtained from PredictDB. We performed a full harmonization of the variants, including recovering strand ambiguous variants. This procedure is discussed in a separate document. (TO-DO: add report that describes harmonization)
LD matrices were computed from a 10% subset of Neale lab subjects. Subjects were matched using the plate and well information from genotyping. We included only biallelic variants with MAF>0.01 in the original Neale Lab GWAS. (TO-DO: add more details [number of subjects, variants, etc])
TO-DO: add enhanced QC reporting (total number of weights, why each variant was missing for all genes)
qclist_all <- list()
qc_files <- paste0(results_dir, "/", list.files(results_dir, pattern="exprqc.Rd"))
for (i in 1:length(qc_files)){
load(qc_files[i])
chr <- unlist(strsplit(rev(unlist(strsplit(qc_files[i], "_")))[1], "[.]"))[1]
qclist_all[[chr]] <- cbind(do.call(rbind, lapply(qclist,unlist)), as.numeric(substring(chr,4)))
}
qclist_all <- data.frame(do.call(rbind, qclist_all))
colnames(qclist_all)[ncol(qclist_all)] <- "chr"
rm(qclist, wgtlist, z_gene_chr)
#number of imputed weights
nrow(qclist_all)
[1] 10896
#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
1070 763 652 417 494 611 548 408 405 434 634 629 195 365 354
16 17 18 19 20 21 22
526 663 160 859 306 114 289
#proportion of imputed weights without missing variants
mean(qclist_all$nmiss==0)
[1] 0.8367291
#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$mu/sample_size #check PVE calculation
#separate gene and SNP results
ctwas_gene_res <- ctwas_res[ctwas_res$type == "gene", ]
ctwas_gene_res <- data.frame(ctwas_gene_res)
ctwas_snp_res <- ctwas_res[ctwas_res$type == "SNP", ]
ctwas_snp_res <- data.frame(ctwas_snp_res)
#add gene information to results
sqlite <- RSQLite::dbDriver("SQLite")
db = RSQLite::dbConnect(sqlite, paste0("/project2/compbio/predictdb/mashr_models/mashr_", weight, ".db"))
query <- function(...) RSQLite::dbGetQuery(db, ...)
gene_info <- query("select gene, genename from extra")
gene_info <- query("select gene, genename, gene_type from extra")
RSQLite::dbDisconnect(db)
ctwas_gene_res <- cbind(ctwas_gene_res, gene_info[sapply(ctwas_gene_res$id, match, gene_info$gene), c("genename", "gene_type")])
#add z score to results
load(paste0(results_dir, "/", analysis_id, "_expr_z_gene.Rd"))
ctwas_gene_res$z <- z_gene[ctwas_gene_res$id,]$z
#load(paste0(results_dir, "/", analysis_id, "_expr_z_snp.Rd")) #for new version, stored after harmonization
z_snp <- readRDS(paste0(results_dir, "/", trait_id, ".RDS")) #for old version, unharmonized
z_snp <- z_snp[z_snp$id %in% ctwas_snp_res$id,] #subset snps to those included in analysis, note some are duplicated, need to match which allele was used
ctwas_snp_res$z <- z_snp$z[match(ctwas_snp_res$id, z_snp$id)] #for duplicated snps, this arbitrarily uses the first allele
ctwas_snp_res$z_flag <- as.numeric(ctwas_snp_res$id %in% z_snp$id[duplicated(z_snp$id)]) #mark the unclear z scores, flag=1
#formatting and rounding for tables
ctwas_gene_res$z <- round(ctwas_gene_res$z,2)
ctwas_snp_res$z <- round(ctwas_snp_res$z,2)
ctwas_gene_res$susie_pip <- round(ctwas_gene_res$susie_pip,3)
ctwas_snp_res$susie_pip <- round(ctwas_snp_res$susie_pip,3)
ctwas_gene_res$mu2 <- round(ctwas_gene_res$mu2,2)
ctwas_snp_res$mu2 <- round(ctwas_snp_res$mu2,2)
ctwas_gene_res$PVE <- signif(ctwas_gene_res$PVE, 2)
ctwas_snp_res$PVE <- signif(ctwas_snp_res$PVE, 2)
#merge gene and snp results with added information
ctwas_gene_res$z_flag=NA
ctwas_snp_res$genename=NA
ctwas_snp_res$gene_type=NA
ctwas_res <- rbind(ctwas_gene_res,
ctwas_snp_res[,colnames(ctwas_gene_res)])
#store columns to report
report_cols <- colnames(ctwas_gene_res)[!(colnames(ctwas_gene_res) %in% c("type", "region_tag1", "region_tag2", "cs_index", "gene_type", "z_flag", "id", "chrom", "pos"))]
first_cols <- c("genename", "region_tag")
report_cols <- c(first_cols, report_cols[!(report_cols %in% first_cols)])
report_cols_snps <- c("id", report_cols[-1])
#get number of SNPs from s1 results; adjust for thin
ctwas_res_s1 <- data.table::fread(paste0(results_dir, "/", analysis_id, "_ctwas.s1.susieIrss.txt"))
n_snps <- sum(ctwas_res_s1$type=="SNP")/thin
rm(ctwas_res_s1)
library(ggplot2)
library(cowplot)
********************************************************
Note: As of version 1.0.0, cowplot does not change the
default ggplot2 theme anymore. To recover the previous
behavior, execute:
theme_set(theme_cowplot())
********************************************************
load(paste0(results_dir, "/", analysis_id, "_ctwas.s2.susieIrssres.Rd"))
df <- data.frame(niter = rep(1:ncol(group_prior_rec), 2),
value = c(group_prior_rec[1,], group_prior_rec[2,]),
group = rep(c("Gene", "SNP"), each = ncol(group_prior_rec)))
df$group <- as.factor(df$group)
df$value[df$group=="SNP"] <- df$value[df$group=="SNP"]*thin #adjust parameter to account for thin argument
p_pi <- ggplot(df, aes(x=niter, y=value, group=group)) +
geom_line(aes(color=group)) +
geom_point(aes(color=group)) +
xlab("Iteration") + ylab(bquote(pi)) +
ggtitle("Prior mean") +
theme_cowplot()
df <- data.frame(niter = rep(1:ncol(group_prior_var_rec), 2),
value = c(group_prior_var_rec[1,], group_prior_var_rec[2,]),
group = rep(c("Gene", "SNP"), each = ncol(group_prior_var_rec)))
df$group <- as.factor(df$group)
p_sigma2 <- ggplot(df, aes(x=niter, y=value, group=group)) +
geom_line(aes(color=group)) +
geom_point(aes(color=group)) +
xlab("Iteration") + ylab(bquote(sigma^2)) +
ggtitle("Prior variance") +
theme_cowplot()
plot_grid(p_pi, p_sigma2)
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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
3.276241e-02 4.791019e-05
#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
5.907624 26.209767
#report sample size
print(sample_size)
[1] 292933
#report group size
group_size <- c(nrow(ctwas_gene_res), n_snps)
print(group_size)
[1] 10896 8696790
#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.007199253 0.037280498
#compare sum(PIP*mu2/sample_size) with above PVE calculation
c(sum(ctwas_gene_res$PVE),sum(ctwas_snp_res$PVE))
[1] 0.0489526 0.5016989
#distribution of PIPs
hist(ctwas_gene_res$susie_pip, xlim=c(0,1), main="Distribution of Gene PIPs")
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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
3212 CCND2 12_4 0.996 26.57 9.0e-05 5.34
7040 INHBB 2_70 0.992 22.43 7.6e-05 4.81
12467 RP11-219B17.3 15_27 0.985 44.27 1.5e-04 7.18
3562 ACVR1C 2_94 0.981 21.16 7.1e-05 4.62
1320 CWF19L1 10_64 0.981 27.59 9.2e-05 -7.09
11790 CYP2A6 19_28 0.975 21.17 7.0e-05 -4.73
5563 ABCG8 2_27 0.966 31.05 1.0e-04 5.88
2359 ABCC3 17_29 0.960 18.71 6.1e-05 4.38
4269 ITGB4 17_42 0.955 20.21 6.6e-05 -4.91
12687 RP4-781K5.7 1_121 0.938 18.26 5.8e-05 -4.17
5978 ZC3H12C 11_65 0.912 18.40 5.7e-05 -4.19
10495 PRMT6 1_66 0.911 24.80 7.7e-05 5.14
2924 EFHD1 2_136 0.910 26.06 8.1e-05 6.05
1146 DNMT3B 20_19 0.910 17.17 5.3e-05 -3.98
1848 CD276 15_35 0.908 32.47 1.0e-04 6.13
1120 CETP 16_31 0.904 18.25 5.6e-05 -4.03
1231 PABPC4 1_24 0.901 20.87 6.4e-05 4.52
1153 TGDS 13_47 0.885 17.41 5.3e-05 -4.00
6682 CYB5R1 1_102 0.868 17.41 5.2e-05 -3.95
5089 SCAF11 12_29 0.856 18.24 5.3e-05 4.13
6093 CSNK1G3 5_75 0.855 21.30 6.2e-05 4.74
10212 IL27 16_23 0.849 21.58 6.3e-05 -4.76
666 COASY 17_25 0.829 17.68 5.0e-05 -3.97
11669 RP11-452H21.4 11_43 0.816 29.71 8.3e-05 5.78
5544 CNIH4 1_114 0.810 17.41 4.8e-05 -3.78
1494 TXN2 22_14 0.803 16.69 4.6e-05 3.67
7547 LIPC 15_26 0.802 18.54 5.1e-05 3.99
#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 |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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
12683 HCP5B 6_24 0.000 10425.85 2.7e-07 -3.52
10663 TRIM31 6_24 0.000 5437.97 9.8e-11 1.07
4833 FLOT1 6_24 0.000 5168.50 8.4e-11 -1.07
1088 USP40 2_137 0.000 3395.97 0.0e+00 -46.64
10651 ABCF1 6_24 0.000 2380.05 7.7e-09 -3.76
5766 PPP1R18 6_24 0.000 2074.46 1.3e-08 -3.94
4836 NRM 6_24 0.000 1159.31 6.1e-12 -0.40
10667 HLA-G 6_24 0.762 1150.05 3.0e-03 -6.69
624 ZNRD1 6_24 0.000 816.36 5.0e-12 0.19
879 DGKD 2_137 0.000 796.64 0.0e+00 -0.07
10747 SLCO1B7 12_16 0.000 624.38 0.0e+00 12.26
3556 HJURP 2_137 0.000 373.30 0.0e+00 10.96
10661 TRIM10 6_24 0.000 346.69 2.0e-12 -0.49
2584 SLCO1B3 12_16 0.000 287.51 0.0e+00 9.93
10648 C6orf136 6_24 0.000 155.03 8.4e-13 0.11
11136 HCG20 6_24 0.000 150.97 4.5e-12 -1.98
4482 SPX 12_16 0.000 150.77 0.0e+00 4.60
10664 RNF39 6_24 0.000 134.08 1.1e-12 -1.02
36 RECQL 12_16 0.000 129.06 0.0e+00 3.81
823 TACR2 10_46 0.014 107.78 5.3e-06 3.90
#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
10667 HLA-G 6_24 0.762 1150.05 3.0e-03 -6.69
12467 RP11-219B17.3 15_27 0.985 44.27 1.5e-04 7.18
5563 ABCG8 2_27 0.966 31.05 1.0e-04 5.88
1848 CD276 15_35 0.908 32.47 1.0e-04 6.13
537 DGAT2 11_42 0.608 46.46 9.6e-05 -7.51
1320 CWF19L1 10_64 0.981 27.59 9.2e-05 -7.09
3212 CCND2 12_4 0.996 26.57 9.0e-05 5.34
11669 RP11-452H21.4 11_43 0.816 29.71 8.3e-05 5.78
2924 EFHD1 2_136 0.910 26.06 8.1e-05 6.05
2004 TGFB1 19_28 0.731 31.30 7.8e-05 5.64
10495 PRMT6 1_66 0.911 24.80 7.7e-05 5.14
7040 INHBB 2_70 0.992 22.43 7.6e-05 4.81
3562 ACVR1C 2_94 0.981 21.16 7.1e-05 4.62
11790 CYP2A6 19_28 0.975 21.17 7.0e-05 -4.73
4269 ITGB4 17_42 0.955 20.21 6.6e-05 -4.91
8142 CNTROB 17_7 0.685 27.81 6.5e-05 -5.71
1231 PABPC4 1_24 0.901 20.87 6.4e-05 4.52
10212 IL27 16_23 0.849 21.58 6.3e-05 -4.76
6093 CSNK1G3 5_75 0.855 21.30 6.2e-05 4.74
2359 ABCC3 17_29 0.960 18.71 6.1e-05 4.38
#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
1088 USP40 2_137 0.000 3395.97 0.0e+00 -46.64
10747 SLCO1B7 12_16 0.000 624.38 0.0e+00 12.26
3556 HJURP 2_137 0.000 373.30 0.0e+00 10.96
8651 MSL2 3_84 0.091 81.97 2.5e-05 10.28
2584 SLCO1B3 12_16 0.000 287.51 0.0e+00 9.93
2586 GOLT1B 12_16 0.000 62.28 0.0e+00 7.53
537 DGAT2 11_42 0.608 46.46 9.6e-05 -7.51
11290 MAPKAPK5-AS1 12_67 0.127 37.01 1.6e-05 -7.21
12467 RP11-219B17.3 15_27 0.985 44.27 1.5e-04 7.18
2541 ALDH2 12_67 0.107 35.10 1.3e-05 7.10
1320 CWF19L1 10_64 0.981 27.59 9.2e-05 -7.09
2536 SH2B3 12_67 0.055 30.20 5.7e-06 6.80
10667 HLA-G 6_24 0.762 1150.05 3.0e-03 -6.69
2170 AHR 7_17 0.052 30.99 5.5e-06 -6.58
4962 EXOC6 10_59 0.099 42.61 1.4e-05 -6.37
8505 HECTD4 12_67 0.209 32.59 2.3e-05 6.33
1191 ERP29 12_67 0.185 29.97 1.9e-05 6.25
10370 TMEM116 12_67 0.185 29.97 1.9e-05 -6.25
9829 ZKSCAN4 6_22 0.095 21.43 6.9e-06 -6.24
10425 AKR1C4 10_6 0.198 31.54 2.1e-05 6.16
#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 |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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 |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#proportion of significant z scores
mean(abs(ctwas_gene_res$z) > sig_thresh)
[1] 0.006240822
#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
1088 USP40 2_137 0.000 3395.97 0.0e+00 -46.64
10747 SLCO1B7 12_16 0.000 624.38 0.0e+00 12.26
3556 HJURP 2_137 0.000 373.30 0.0e+00 10.96
8651 MSL2 3_84 0.091 81.97 2.5e-05 10.28
2584 SLCO1B3 12_16 0.000 287.51 0.0e+00 9.93
2586 GOLT1B 12_16 0.000 62.28 0.0e+00 7.53
537 DGAT2 11_42 0.608 46.46 9.6e-05 -7.51
11290 MAPKAPK5-AS1 12_67 0.127 37.01 1.6e-05 -7.21
12467 RP11-219B17.3 15_27 0.985 44.27 1.5e-04 7.18
2541 ALDH2 12_67 0.107 35.10 1.3e-05 7.10
1320 CWF19L1 10_64 0.981 27.59 9.2e-05 -7.09
2536 SH2B3 12_67 0.055 30.20 5.7e-06 6.80
10667 HLA-G 6_24 0.762 1150.05 3.0e-03 -6.69
2170 AHR 7_17 0.052 30.99 5.5e-06 -6.58
4962 EXOC6 10_59 0.099 42.61 1.4e-05 -6.37
8505 HECTD4 12_67 0.209 32.59 2.3e-05 6.33
1191 ERP29 12_67 0.185 29.97 1.9e-05 6.25
10370 TMEM116 12_67 0.185 29.97 1.9e-05 -6.25
9829 ZKSCAN4 6_22 0.095 21.43 6.9e-06 -6.24
10425 AKR1C4 10_6 0.198 31.54 2.1e-05 6.16
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: 2_137"
genename region_tag susie_pip mu2 PVE z
10567 GIGYF2 2_137 0 54.05 0 -5.12
9340 C2orf82 2_137 0 5.37 0 0.45
620 NGEF 2_137 0 44.21 0 2.52
8006 INPP5D 2_137 0 36.50 0 4.05
879 DGKD 2_137 0 796.64 0 -0.07
1088 USP40 2_137 0 3395.97 0 -46.64
3556 HJURP 2_137 0 373.30 0 10.96
11098 AC006037.2 2_137 0 15.79 0 1.50
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
[1] "Region: 12_16"
genename region_tag susie_pip mu2 PVE z
2584 SLCO1B3 12_16 0 287.51 0 9.93
10747 SLCO1B7 12_16 0 624.38 0 12.26
3400 IAPP 12_16 0 95.90 0 5.16
3399 PYROXD1 12_16 0 11.90 0 0.89
36 RECQL 12_16 0 129.06 0 3.81
2586 GOLT1B 12_16 0 62.28 0 7.53
4482 SPX 12_16 0 150.77 0 4.60
2587 LDHB 12_16 0 8.40 0 0.07
3401 KCNJ8 12_16 0 6.08 0 -0.33
689 ABCC9 12_16 0 7.00 0 2.11
2590 C2CD5 12_16 0 10.33 0 -1.39
5073 ETNK1 12_16 0 9.07 0 0.30
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
[1] "Region: 3_84"
genename region_tag susie_pip mu2 PVE z
796 PPP2R3A 3_84 0.145 13.81 6.9e-06 -2.46
8651 MSL2 3_84 0.091 81.97 2.5e-05 10.28
2795 PCCB 3_84 0.064 5.48 1.2e-06 1.22
3148 STAG1 3_84 0.066 4.63 1.0e-06 -0.03
6584 NCK1 3_84 0.064 7.94 1.7e-06 -2.19
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
[1] "Region: 11_42"
genename region_tag susie_pip mu2 PVE z
7611 XRRA1 11_42 0.120 7.42 3.0e-06 -1.00
3170 SPCS2 11_42 0.095 5.54 1.8e-06 0.73
6901 NEU3 11_42 0.085 4.45 1.3e-06 0.39
4848 SLCO2B1 11_42 0.085 4.44 1.3e-06 -0.22
12001 TPBGL 11_42 0.110 6.57 2.5e-06 0.68
6617 GDPD5 11_42 0.120 9.34 3.8e-06 1.95
8328 MAP6 11_42 0.098 5.57 1.9e-06 0.24
7603 MOGAT2 11_42 0.264 12.77 1.2e-05 0.85
537 DGAT2 11_42 0.608 46.46 9.6e-05 -7.51
10381 UVRAG 11_42 0.095 6.40 2.1e-06 1.62
1082 WNT11 11_42 0.117 7.35 2.9e-06 1.06
11773 RP11-619A14.3 11_42 0.095 5.41 1.8e-06 0.61
4849 THAP12 11_42 0.091 5.17 1.6e-06 -0.66
12265 RP11-111M22.5 11_42 0.099 5.88 2.0e-06 0.80
11766 RP11-111M22.3 11_42 0.085 4.48 1.3e-06 0.31
11751 RP11-672A2.4 11_42 0.092 5.08 1.6e-06 0.53
9350 TSKU 11_42 0.089 4.75 1.4e-06 0.25
905 ACER3 11_42 0.265 14.43 1.3e-05 -2.04
5976 CAPN5 11_42 0.196 11.97 8.0e-06 1.72
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
[1] "Region: 12_67"
genename region_tag susie_pip mu2 PVE z
5112 TCHP 12_67 0.495 15.40 2.6e-05 3.21
5111 GIT2 12_67 0.297 14.18 1.4e-05 3.36
8639 C12orf76 12_67 0.069 6.01 1.4e-06 -0.30
3515 IFT81 12_67 0.074 9.67 2.4e-06 2.50
10093 ANAPC7 12_67 0.072 7.73 1.9e-06 2.16
2531 ARPC3 12_67 0.078 6.76 1.8e-06 0.54
10684 FAM216A 12_67 0.063 7.90 1.7e-06 2.42
2532 GPN3 12_67 0.055 4.79 8.9e-07 1.40
2533 VPS29 12_67 0.055 4.82 9.0e-07 -1.41
10683 TCTN1 12_67 0.069 5.68 1.3e-06 0.02
3517 HVCN1 12_67 0.108 11.91 4.4e-06 2.67
9717 PPP1CC 12_67 0.108 11.85 4.4e-06 -2.65
10375 FAM109A 12_67 0.057 6.81 1.3e-06 -1.46
2536 SH2B3 12_67 0.055 30.20 5.7e-06 6.80
10680 ATXN2 12_67 0.057 17.13 3.3e-06 3.97
2541 ALDH2 12_67 0.107 35.10 1.3e-05 7.10
11290 MAPKAPK5-AS1 12_67 0.127 37.01 1.6e-05 -7.21
1191 ERP29 12_67 0.185 29.97 1.9e-05 6.25
10370 TMEM116 12_67 0.185 29.97 1.9e-05 -6.25
2544 NAA25 12_67 0.155 27.83 1.5e-05 -6.12
8505 HECTD4 12_67 0.209 32.59 2.3e-05 6.33
9084 PTPN11 12_67 0.072 6.64 1.6e-06 -0.78
Version | Author | Date |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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
326005 rs72834643 6_20 1.000 118.88 4.1e-04 9.73
326026 rs115740542 6_20 1.000 180.46 6.2e-04 12.66
370129 rs12208357 6_103 1.000 80.73 2.8e-04 -6.65
381915 rs542176135 7_17 1.000 152.94 5.2e-04 -8.38
533135 rs6480402 10_46 1.000 469.87 1.6e-03 8.90
603061 rs11045819 12_16 1.000 2552.27 8.7e-03 -14.34
603078 rs4363657 12_16 1.000 1590.67 5.4e-03 43.78
797812 rs59616136 19_14 1.000 49.81 1.7e-04 7.00
897096 rs13031505 2_136 1.000 60.72 2.1e-04 -8.07
900079 rs13402456 2_137 1.000 3952.27 1.3e-02 -44.40
900100 rs28948388 2_137 1.000 2663.80 9.1e-03 -35.74
900242 rs6431631 2_137 1.000 8372.22 2.9e-02 49.75
900289 rs150783152 2_137 1.000 3828.43 1.3e-02 45.10
909053 rs1611236 6_24 1.000 27586.62 9.4e-02 -3.60
326868 rs3130253 6_23 0.999 41.11 1.4e-04 -6.63
508420 rs115478735 9_70 0.999 65.26 2.2e-04 -8.02
381937 rs4721597 7_17 0.996 89.08 3.0e-04 1.94
554431 rs76153913 11_2 0.996 47.07 1.6e-04 6.70
805958 rs113345881 19_32 0.996 35.92 1.2e-04 6.12
629891 rs653178 12_67 0.995 155.80 5.3e-04 -13.12
897033 rs34247311 2_136 0.995 38.53 1.3e-04 4.44
613305 rs7397189 12_36 0.993 35.11 1.2e-04 5.69
817689 rs34507316 20_13 0.993 33.09 1.1e-04 5.60
805956 rs814573 19_32 0.989 39.63 1.3e-04 -6.68
36470 rs2779116 1_78 0.981 68.79 2.3e-04 -8.25
744704 rs2608604 16_53 0.976 55.78 1.9e-04 -6.33
796041 rs141645070 19_10 0.969 30.48 1.0e-04 -5.27
798597 rs3794991 19_15 0.952 131.75 4.3e-04 11.64
805892 rs1551891 19_32 0.946 35.90 1.2e-04 7.88
845995 rs34662558 22_10 0.935 30.92 9.9e-05 -5.19
602776 rs7962574 12_15 0.922 45.95 1.4e-04 -8.40
805961 rs12721109 19_32 0.918 33.02 1.0e-04 6.44
494338 rs9410381 9_45 0.913 74.54 2.3e-04 8.62
233289 rs17238095 4_72 0.910 29.26 9.1e-05 5.18
436846 rs12549737 8_24 0.896 30.05 9.2e-05 5.15
508805 rs34755157 9_71 0.896 29.05 8.9e-05 -5.03
900250 rs35025179 2_137 0.890 5886.79 1.8e-02 -49.78
558884 rs34623292 11_10 0.888 30.00 9.1e-05 -5.06
602781 rs73080739 12_15 0.879 34.16 1.0e-04 -7.24
325844 rs75080831 6_19 0.854 60.06 1.8e-04 7.96
602816 rs10770693 12_15 0.851 59.30 1.7e-04 8.86
294604 rs4566840 5_66 0.840 31.47 9.0e-05 -5.46
557453 rs4910498 11_8 0.833 50.56 1.4e-04 6.71
#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 |
---|---|---|
b14741c | wesleycrouse | 2021-09-06 |
#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
909053 rs1611236 6_24 1.000 27586.62 9.4e-02 -3.60
909076 rs1611248 6_24 0.237 27585.63 2.2e-02 -3.63
909080 rs1611252 6_24 0.140 27585.51 1.3e-02 -3.62
909035 rs111734624 6_24 0.146 27585.47 1.4e-02 -3.62
909032 rs2508055 6_24 0.146 27585.46 1.4e-02 -3.62
909097 rs1611260 6_24 0.138 27585.44 1.3e-02 -3.62
909103 rs1611265 6_24 0.137 27585.42 1.3e-02 -3.62
908971 rs1633033 6_24 0.253 27585.30 2.4e-02 -3.63
908984 rs2844838 6_24 0.129 27584.82 1.2e-02 -3.63
909105 rs1611267 6_24 0.057 27584.79 5.4e-03 -3.62
909107 rs2394171 6_24 0.049 27584.76 4.6e-03 -3.61
909109 rs2893981 6_24 0.049 27584.71 4.6e-03 -3.61
909039 rs1611228 6_24 0.063 27584.64 5.9e-03 -3.62
909028 rs1737020 6_24 0.020 27584.20 1.8e-03 -3.60
909029 rs1737019 6_24 0.020 27584.20 1.8e-03 -3.60
908988 rs1633032 6_24 0.066 27582.93 6.2e-03 -3.62
909022 rs1633020 6_24 0.026 27581.16 2.4e-03 -3.61
909026 rs1633018 6_24 0.009 27580.58 8.2e-04 -3.60
909051 rs1611234 6_24 0.009 27578.97 8.1e-04 -3.60
908911 rs1610726 6_24 0.014 27577.73 1.3e-03 -3.61
908979 rs2844840 6_24 0.034 27575.49 3.2e-03 -3.62
909306 rs3129185 6_24 0.041 27575.16 3.9e-03 -3.63
909321 rs1736999 6_24 0.029 27573.96 2.8e-03 -3.63
909334 rs1633001 6_24 0.055 27572.06 5.2e-03 -3.64
909074 rs1611246 6_24 0.014 27570.41 1.4e-03 -3.62
909510 rs1632980 6_24 0.025 27570.15 2.4e-03 -3.63
909007 rs1614309 6_24 0.016 27562.78 1.5e-03 -3.63
909006 rs1633030 6_24 0.000 27539.27 5.5e-06 -3.61
909119 rs9258382 6_24 0.001 27512.70 9.1e-05 -3.69
909116 rs9258379 6_24 0.000 27465.77 3.9e-09 -3.62
909065 rs1611241 6_24 0.004 27439.03 4.0e-04 -3.84
909010 rs1633028 6_24 0.000 27397.93 5.6e-09 -3.72
909068 rs1611244 6_24 0.000 27291.46 3.4e-09 -3.64
909023 rs2735042 6_24 0.000 27248.96 3.5e-09 -3.61
909104 rs1611266 6_24 0.000 27051.57 1.6e-08 -3.97
909077 rs1611249 6_24 0.000 26937.56 3.0e-08 -4.10
909043 rs1611230 6_24 0.000 26870.84 2.6e-08 -4.08
909092 rs145043018 6_24 0.000 26865.68 2.8e-08 -4.09
909102 rs147376303 6_24 0.000 26865.55 2.8e-08 -4.09
909113 rs9258376 6_24 0.000 26865.21 2.8e-08 -4.09
909120 rs1633016 6_24 0.000 26864.88 2.8e-08 -4.09
908965 rs1633035 6_24 0.000 26864.28 2.8e-08 -4.09
909173 rs1633014 6_24 0.000 26862.45 2.7e-08 -4.09
908998 rs1618670 6_24 0.000 26862.21 2.7e-08 -4.09
909025 rs1633019 6_24 0.000 26859.77 2.4e-08 -4.06
909287 rs1633010 6_24 0.000 26854.25 2.7e-08 -4.09
909411 rs909722 6_24 0.000 26849.65 2.6e-08 -4.08
909443 rs1610713 6_24 0.000 26849.60 2.6e-08 -4.09
909368 rs1736991 6_24 0.000 26849.44 2.7e-08 -4.09
909423 rs1610648 6_24 0.000 26848.15 2.8e-08 -4.10
#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
909053 rs1611236 6_24 1.000 27586.62 0.09400 -3.60
900242 rs6431631 2_137 1.000 8372.22 0.02900 49.75
908971 rs1633033 6_24 0.253 27585.30 0.02400 -3.63
909076 rs1611248 6_24 0.237 27585.63 0.02200 -3.63
900250 rs35025179 2_137 0.890 5886.79 0.01800 -49.78
909032 rs2508055 6_24 0.146 27585.46 0.01400 -3.62
909035 rs111734624 6_24 0.146 27585.47 0.01400 -3.62
900079 rs13402456 2_137 1.000 3952.27 0.01300 -44.40
900289 rs150783152 2_137 1.000 3828.43 0.01300 45.10
909080 rs1611252 6_24 0.140 27585.51 0.01300 -3.62
909097 rs1611260 6_24 0.138 27585.44 0.01300 -3.62
909103 rs1611265 6_24 0.137 27585.42 0.01300 -3.62
908984 rs2844838 6_24 0.129 27584.82 0.01200 -3.63
900100 rs28948388 2_137 1.000 2663.80 0.00910 -35.74
603061 rs11045819 12_16 1.000 2552.27 0.00870 -14.34
908988 rs1633032 6_24 0.066 27582.93 0.00620 -3.62
602962 rs12366506 12_16 0.714 2465.44 0.00600 23.44
909039 rs1611228 6_24 0.063 27584.64 0.00590 -3.62
603078 rs4363657 12_16 1.000 1590.67 0.00540 43.78
909105 rs1611267 6_24 0.057 27584.79 0.00540 -3.62
909334 rs1633001 6_24 0.055 27572.06 0.00520 -3.64
909107 rs2394171 6_24 0.049 27584.76 0.00460 -3.61
909109 rs2893981 6_24 0.049 27584.71 0.00460 -3.61
602968 rs11045612 12_16 0.532 2463.24 0.00450 23.44
909306 rs3129185 6_24 0.041 27575.16 0.00390 -3.63
908979 rs2844840 6_24 0.034 27575.49 0.00320 -3.62
909321 rs1736999 6_24 0.029 27573.96 0.00280 -3.63
909022 rs1633020 6_24 0.026 27581.16 0.00240 -3.61
909510 rs1632980 6_24 0.025 27570.15 0.00240 -3.63
900253 rs17863805 2_137 0.103 5879.66 0.00210 -49.75
602979 rs73062442 12_16 0.215 2463.49 0.00180 23.39
909028 rs1737020 6_24 0.020 27584.20 0.00180 -3.60
909029 rs1737019 6_24 0.020 27584.20 0.00180 -3.60
602971 rs11513221 12_16 0.198 2463.53 0.00170 23.38
533135 rs6480402 10_46 1.000 469.87 0.00160 8.90
909007 rs1614309 6_24 0.016 27562.78 0.00150 -3.63
909074 rs1611246 6_24 0.014 27570.41 0.00140 -3.62
908911 rs1610726 6_24 0.014 27577.73 0.00130 -3.61
909026 rs1633018 6_24 0.009 27580.58 0.00082 -3.60
909051 rs1611234 6_24 0.009 27578.97 0.00081 -3.60
533140 rs35233497 10_46 0.525 352.95 0.00063 -4.24
326026 rs115740542 6_20 1.000 180.46 0.00062 12.66
533143 rs79086908 10_46 0.474 352.88 0.00057 -4.24
629891 rs653178 12_67 0.995 155.80 0.00053 -13.12
381915 rs542176135 7_17 1.000 152.94 0.00052 -8.38
602872 rs3060461 12_16 0.317 477.25 0.00052 -21.86
914179 rs3132129 6_24 0.094 1434.35 0.00046 -7.09
602876 rs35290079 12_16 0.279 474.92 0.00045 -21.76
798597 rs3794991 19_15 0.952 131.75 0.00043 11.64
326005 rs72834643 6_20 1.000 118.88 0.00041 9.73
#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
900250 rs35025179 2_137 0.890 5886.79 1.8e-02 -49.78
900242 rs6431631 2_137 1.000 8372.22 2.9e-02 49.75
900253 rs17863805 2_137 0.103 5879.66 2.1e-03 -49.75
900246 rs9287649 2_137 0.001 5867.77 1.0e-05 -49.74
900245 rs4663973 2_137 0.006 5873.92 1.3e-04 -49.71
900244 rs4663335 2_137 0.000 5862.95 1.3e-07 -49.63
899151 rs62192912 2_137 0.000 3896.63 0.0e+00 -49.52
900224 rs1500482 2_137 0.000 8260.81 0.0e+00 49.20
899798 rs72980395 2_137 0.000 6695.23 0.0e+00 48.79
900267 rs12468017 2_137 0.000 5622.20 0.0e+00 -48.71
900247 rs6431632 2_137 0.000 8103.88 0.0e+00 48.44
900260 rs17862885 2_137 0.000 5661.36 0.0e+00 -48.44
900220 rs8330 2_137 0.000 8057.30 0.0e+00 48.40
900249 rs10199525 2_137 0.000 8103.55 0.0e+00 48.32
899649 rs7584554 2_137 0.000 4661.10 0.0e+00 48.31
900262 rs116149023 2_137 0.000 5614.64 0.0e+00 -48.31
900254 rs11903524 2_137 0.000 8101.74 0.0e+00 48.30
900248 rs6431633 2_137 0.000 8100.21 0.0e+00 48.29
900251 rs7583811 2_137 0.000 8100.57 0.0e+00 48.29
900257 rs7596785 2_137 0.000 8097.59 0.0e+00 48.26
899999 rs62192806 2_137 0.000 6434.77 0.0e+00 48.16
900258 rs10929304 2_137 0.000 8076.73 0.0e+00 48.13
899620 rs838716 2_137 0.000 2530.54 0.0e+00 48.10
899968 rs62192778 2_137 0.000 4587.29 0.0e+00 -48.02
899602 rs838711 2_137 0.000 2532.49 0.0e+00 47.95
899956 rs62192777 2_137 0.000 4568.29 0.0e+00 -47.92
899596 rs838709 2_137 0.000 2502.54 0.0e+00 47.89
899595 rs838708 2_137 0.000 2499.08 0.0e+00 47.88
900329 rs6716168 2_137 0.000 7644.44 0.0e+00 -47.87
900221 rs17862880 2_137 0.000 5324.94 0.0e+00 -47.84
900206 rs55891750 2_137 0.000 10773.35 0.0e+00 47.71
900326 rs730673 2_137 0.000 7554.83 0.0e+00 -47.70
899842 rs142848041 2_137 0.000 6644.41 0.0e+00 47.35
900264 rs80334668 2_137 0.000 6138.56 0.0e+00 47.35
900226 rs11563250 2_137 0.000 5271.27 0.0e+00 -47.28
900227 rs28900409 2_137 0.000 5270.77 0.0e+00 -47.26
899819 rs139170694 2_137 0.000 6586.89 0.0e+00 47.16
900109 rs12469671 2_137 0.000 4785.04 0.0e+00 -47.07
899793 rs144648801 2_137 0.000 7414.14 0.0e+00 47.00
900255 rs561956448 2_137 0.000 5211.67 0.0e+00 -47.00
900218 rs10929303 2_137 0.000 7618.29 0.0e+00 46.69
899559 rs6727030 2_137 0.000 2362.38 0.0e+00 46.64
899732 rs76218013 2_137 0.000 5132.55 0.0e+00 46.51
899544 rs1550532 2_137 0.000 2338.27 0.0e+00 46.49
899251 rs62192932 2_137 0.000 3449.63 0.0e+00 46.41
900231 rs7586006 2_137 0.000 5197.26 0.0e+00 -46.39
900229 rs10199882 2_137 0.000 5200.72 0.0e+00 -46.38
899704 rs11689634 2_137 0.000 7211.50 0.0e+00 46.34
899738 rs72980357 2_137 0.000 7107.81 0.0e+00 46.06
899983 rs78198652 2_137 0.000 3468.90 0.0e+00 -45.92
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] cowplot_1.0.0 ggplot2_3.3.3
loaded via a namespace (and not attached):
[1] Biobase_2.44.0 httr_1.4.1
[3] bit64_4.0.5 assertthat_0.2.1
[5] stats4_3.6.1 blob_1.2.1
[7] BSgenome_1.52.0 GenomeInfoDbData_1.2.1
[9] Rsamtools_2.0.0 yaml_2.2.0
[11] progress_1.2.2 pillar_1.6.1
[13] RSQLite_2.2.7 lattice_0.20-38
[15] glue_1.4.2 digest_0.6.20
[17] GenomicRanges_1.36.0 promises_1.0.1
[19] XVector_0.24.0 colorspace_1.4-1
[21] htmltools_0.3.6 httpuv_1.5.1
[23] Matrix_1.2-18 XML_3.98-1.20
[25] pkgconfig_2.0.3 biomaRt_2.40.1
[27] zlibbioc_1.30.0 purrr_0.3.4
[29] scales_1.1.0 whisker_0.3-2
[31] later_0.8.0 BiocParallel_1.18.0
[33] git2r_0.26.1 tibble_3.1.2
[35] farver_2.1.0 generics_0.0.2
[37] IRanges_2.18.1 ellipsis_0.3.2
[39] withr_2.4.1 cachem_1.0.5
[41] SummarizedExperiment_1.14.1 GenomicFeatures_1.36.3
[43] BiocGenerics_0.30.0 magrittr_2.0.1
[45] crayon_1.4.1 memoise_2.0.0
[47] evaluate_0.14 fs_1.3.1
[49] fansi_0.5.0 tools_3.6.1
[51] data.table_1.14.0 prettyunits_1.0.2
[53] hms_1.1.0 lifecycle_1.0.0
[55] matrixStats_0.57.0 stringr_1.4.0
[57] S4Vectors_0.22.1 munsell_0.5.0
[59] DelayedArray_0.10.0 AnnotationDbi_1.46.0
[61] Biostrings_2.52.0 compiler_3.6.1
[63] GenomeInfoDb_1.20.0 rlang_0.4.11
[65] grid_3.6.1 RCurl_1.98-1.1
[67] VariantAnnotation_1.30.1 labeling_0.3
[69] bitops_1.0-6 rmarkdown_1.13
[71] gtable_0.3.0 DBI_1.1.1
[73] R6_2.5.0 GenomicAlignments_1.20.1
[75] dplyr_1.0.7 knitr_1.23
[77] rtracklayer_1.44.0 utf8_1.2.1
[79] fastmap_1.1.0 bit_4.0.4
[81] workflowr_1.6.2 rprojroot_2.0.2
[83] stringi_1.4.3 parallel_3.6.1
[85] Rcpp_1.0.6 vctrs_0.3.8
[87] tidyselect_1.1.0 xfun_0.8