Progress report 3-28-2016
Figure 1. Weight effect of rho
Figure 2. Weight effect of lung P value
library(ggplot2)## Warning: package 'ggplot2' was built under R version 3.2.4library(reshape)
library(knitr)
library(xtable)merged.mouse.eQTL.min.variance2<-read.table(file="2016-02-22-merged.mouse.eQTL.min.variance.txt", header=T)
merged.mouse.eQTL.min.variance2$abs_liver.beta<-abs(merged.mouse.eQTL.min.variance2$liver.beta)
merged.mouse.eQTL.min.variance2$abs_lung.beta<-abs(merged.mouse.eQTL.min.variance2$lung.beta)
merged.mouse.eQTL.min.variance2$neg_log_lung_pvalue<--log(merged.mouse.eQTL.min.variance2$lung_pvalue)
merged.mouse.eQTL<-merged.mouse.eQTL.min.variance2
markers<-merged.mouse.eQTL[, 1]
# Yg=Ag + Bg*Xsnp+V
betas.hat<-abs(merged.mouse.eQTL[,2])
se<-(merged.mouse.eQTL[,3])
liver_residual_variance<-merged.mouse.eQTL[,4]
Z<-as.matrix(merged.mouse.eQTL[,10])
Z<-replace(Z,is.na(Z),0)
Z<-data.frame(1,Z)
Z<-as.matrix(Z)
rowLength<-length(markers)
lmsummary<-summary(lm(abs_liver.beta~-1+Z, data=merged.mouse.eQTL))
tau<-lmsummary$sigma**2
tau## [1] 0.009699749gamma<-as.matrix(lmsummary$coefficients[,1])
Z_transpose<-t(Z)
identity<-diag(nrow=rowLength)
betas.hat<-as.matrix(betas.hat)
V <- matrix(0, rowLength, rowLength)
diag(V) <-liver_residual_variance
Tau<- diag(tau, rowLength, rowLength)
s<-V + Taudiag.inverse <- function(x){diag(1/diag(x), nrow(x), ncol(x))}
diag.multi <- function(x,y){diag(diag(x)*diag(y), nrow(x), ncol(x))}
S <-diag.inverse(s)
omega<-diag.multi(S, V)
omega.diag<-diag(omega )
# betas.thea<- S %*% Z %*% gamma + (identity-S) %*% betas.hat
betas.tieda<- omega %*% Z %*% gamma + (identity-omega) %*% betas.hat
# crbetas.tieda<- cromega %*% Z %*% gamma + (identity-cromega) %*% betas.hat
regbeta <-Z %*% gammaliver.mouse.eQTL.bayesian.4tau <- read.table(file="2016-03-15_liver.mouse.eQTL.bayesian.4tau.txt", head=TRUE)liver.mouse.eQTL.bayesian.4tau$fzm <- -log10(exp(-liver.mouse.eQTL.bayesian.4tau$neg_log_lung_pvalue))
range(liver.mouse.eQTL.bayesian.4tau$fzm)## [1] 8.678947e-16 1.727677e+01liver.mouse.eQTL.bayesian.4tau$ratio_fzm <- max(liver.mouse.eQTL.bayesian.4tau$fzm)/liver.mouse.eQTL.bayesian.4tau$fzm
range(liver.mouse.eQTL.bayesian.4tau$ratio_fzm)## [1] 1.000000e+00 1.990653e+16rho.optimization <- matrix(0, nrow=nrow(liver.mouse.eQTL.bayesian.4tau), ncol=7)
colnames(rho.optimization)<-c("rho","tmm","tau", "omega","beta_tieda", "n.betas.tieda.se","p.below.0" )
nomega.diag<-diag(nomega )
rho <- seq(1,51, by=10)*tau
rho
result <- NULL
for (i in 1:length(rho)) {
rho.optimization[ ,1] <- rho[i]
rho.optimization[ ,2] <- (rho[i]/tau)^liver.mouse.eQTL.bayesian.4tau$ratio_fzm
rho.optimization[ ,3] <-tau*((rho[i]/tau)^liver.mouse.eQTL.bayesian.4tau$ratio_fzm)
nTau<- diag(rho.optimization[i,3], rowLength, rowLength)
ns<-V + nTau
nS <- diag.inverse(ns)
nomega<-diag.multi(nS, V)
rho.optimization[ ,4] <- diag(nomega )
rho.optimization[ ,5] <-nomega %*% Z %*% gamma + (identity-nomega) %*% betas.hat
nTau_invert<-diag.inverse(nTau)
V_invert<-diag.inverse(V)
nPS_invert<-nTau_invert+ diag.multi(diag.multi(V_invert, Z), Z_transpose)
nPS<-diag.inverse(nPS_invert)
nps<-diag(nPS)
nps.long <- melt(nps)
rho.optimization[ ,6] <-(nps.long$value)^0.5
rho.optimization[ ,7] <- pnorm(0, rho.optimization[ ,5], rho.optimization[ ,6])
result <- rbind(result,rho.optimization)
}
result$p.below.0 <- pnorm(0, rho.optimization[ ,5], rho.optimization[ ,6])
write.table(result, file="2016-03-22_liver.mouse.eQTL.bayesian.result.txt",col.names=TRUE,row.names=FALSE,quote=FALSE)liver.mouse.eQTL.bayesian.result <- read.table(file="2016-03-22_liver.mouse.eQTL.bayesian.result.txt", header=T)
knitr::kable(head(liver.mouse.eQTL.bayesian.result))| rho | tmm | tau | omega | beta_tieda | n.betas.tieda.se | p.below.0 |
|---|---|---|---|---|---|---|
| 0.0096997 | 1 | 0.0096997 | 0.6609725 | 0.0335649 | 0.0800704 | 0.3375376 |
| 0.0096997 | 1 | 0.0096997 | 0.2791983 | 0.0236274 | 0.0983765 | 0.4050981 |
| 0.0096997 | 1 | 0.0096997 | 0.3255088 | 0.0140638 | 0.0561903 | 0.4011821 |
| 0.0096997 | 1 | 0.0096997 | 0.1586587 | 0.0416051 | 0.0982601 | 0.3359952 |
| 0.0096997 | 1 | 0.0096997 | 0.8895190 | 0.0401108 | 0.0928876 | 0.3329360 |
| 0.0096997 | 1 | 0.0096997 | 0.3335208 | 0.0306453 | 0.0984015 | 0.3777363 |
knitr::kable(tail(liver.mouse.eQTL.bayesian.result))| rho | tmm | tau | omega | beta_tieda | n.betas.tieda.se | p.below.0 | |
|---|---|---|---|---|---|---|---|
| 66049 | 0.4946872 | 1.180091 | 0.0114466 | 0.6634209 | 0.0473521 | 0.1046325 | 0.3254342 |
| 66050 | 0.4946872 | 1.472411 | 0.0142820 | 0.3777155 | 0.0383391 | 0.0643198 | 0.2755646 |
| 66051 | 0.4946872 | 1.096974 | 0.0106404 | 0.3395461 | 0.0197970 | 0.1045669 | 0.4249194 |
| 66052 | 0.4946872 | 1.139850 | 0.0110563 | 0.5864554 | 0.0408051 | 0.0801457 | 0.3053284 |
| 66053 | 0.4946872 | 1.306410 | 0.0126718 | 0.4981594 | 0.0320585 | 0.1046096 | 0.3796280 |
| 66054 | 0.4946872 | 1.296591 | 0.0125766 | 0.6772621 | 0.0399182 | 0.0861273 | 0.3215104 |
result.df <-liver.mouse.eQTL.bayesian.result
result.df$rho.class <- factor(result.df$rho/tau)a <-liver.mouse.eQTL.bayesian.4tau[, c(1:2, 6, 7, 9)]
a <-rbind(a, a, a, a, a,a)
dim(a)## [1] 66054 5new.result.df<-cbind(a, result.df)knitr::kable(head(liver.mouse.eQTL.bayesian.4tau)) | ensembl_id | betas.hat | betas.hat.se | betas.tieda | betas.tieda.se | liver_residual_variance | liver_pvalue | abs_lung.beta | neg_log_lung_pvalue | p.below.0 | betas.tieda2m | betas.tieda3rd | betas.tieda4max | fzm | ratio_fzm |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| ENSMUSG00000000001 | 0.0159861 | 0.0256247 | 0.0335649 | 0.0800559 | 0.0189108 | 0.5377718 | 0.0326900 | 2.0960974 | 0.3375098 | 0.0328346 | 0.0277012 | 0.0161111 | 0.9103235 | 18.97871 |
| ENSMUSG00000000003 | 0.0168807 | 0.0126533 | 0.0236274 | 0.0520236 | 0.0037571 | 0.1929225 | 0.0295496 | 2.7653413 | 0.3248544 | 0.0230580 | 0.0201492 | 0.0169033 | 1.2009724 | 14.38565 |
| ENSMUSG00000000037 | 0.0038651 | 0.0136293 | 0.0140638 | 0.0561845 | 0.0046811 | 0.7788138 | 0.0175948 | 1.1814994 | 0.4011720 | 0.0132541 | 0.0089752 | 0.0039017 | 0.5131187 | 33.67012 |
| ENSMUSG00000000049 | 0.0422333 | 0.0078085 | 0.0416051 | 0.0392200 | 0.0018292 | 0.0000091 | 0.0238839 | 0.2321777 | 0.1443878 | 0.0416661 | 0.0419531 | 0.0422315 | 0.1008335 | 171.33960 |
| ENSMUSG00000000056 | 0.0636041 | 0.0514812 | 0.0401108 | 0.0928852 | 0.0780959 | 0.2269175 | 0.0216750 | 1.5218109 | 0.3329319 | 0.0404380 | 0.0434054 | 0.0630989 | 0.6609141 | 26.14072 |
| ENSMUSG00000000058 | 0.0256111 | 0.0129823 | 0.0306453 | 0.0568619 | 0.0048540 | 0.0584692 | 0.0288547 | 1.2292565 | 0.2949631 | 0.0302500 | 0.0281486 | 0.0256294 | 0.5338593 | 32.36202 |
new.result.df2 <- new.result.df
by(new.result.df2, new.result.df2[, "rho.class"], head)## new.result.df2[, "rho.class"]: 1
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 1 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 2 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 3 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 4 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 5 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 6 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega beta_tieda
## 1 2.0960974 0.009699749 1 0.009699749 0.6609725 0.03356494
## 2 2.7653413 0.009699749 1 0.009699749 0.2791983 0.02362737
## 3 1.1814994 0.009699749 1 0.009699749 0.3255088 0.01406380
## 4 0.2321777 0.009699749 1 0.009699749 0.1586587 0.04160507
## 5 1.5218109 0.009699749 1 0.009699749 0.8895190 0.04011077
## 6 1.2292565 0.009699749 1 0.009699749 0.3335208 0.03064532
## n.betas.tieda.se p.below.0 rho.class
## 1 0.08007039 0.3375376 1
## 2 0.09837648 0.4050981 1
## 3 0.05619033 0.4011821 1
## 4 0.09826008 0.3359952 1
## 5 0.09288763 0.3329360 1
## 6 0.09840149 0.3777363 1
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 11
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 11010 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 11011 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 11012 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 11013 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 11014 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 11015 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega
## 11010 2.0960974 0.1066972 1.134675 0.011006067 0.6226813
## 11011 2.7653413 0.1066972 1.181384 0.011459129 0.2469161
## 11012 1.1814994 0.1066972 1.073815 0.010415732 0.2900264
## 11013 0.2321777 0.1066972 1.014093 0.009836452 0.1376522
## 11014 1.5218109 0.1066972 1.096069 0.010631596 0.8720437
## 11015 1.2292565 0.1066972 1.076910 0.010445758 0.2975504
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 11010 0.03254657 0.08447121 0.3500084 11
## 11011 0.02284729 0.10690506 0.4153844 11
## 11012 0.01295208 0.05764937 0.4111180 11
## 11013 0.04168825 0.10675573 0.3480830 11
## 11014 0.04057231 0.09996430 0.3424198 11
## 11015 0.03010237 0.10693716 0.3891650 11
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 21
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 22019 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 22020 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 22021 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 22022 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 22023 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 22024 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega
## 22019 2.0960974 0.2036947 1.174001 0.011387518 0.6404250
## 22020 2.7653413 0.2036947 1.235698 0.011985962 0.2613693
## 22021 1.1814994 0.2036947 1.094636 0.010617697 0.3059778
## 22022 0.2321777 0.2036947 1.017928 0.009873644 0.1469577
## 22023 1.5218109 0.2036947 1.123520 0.010897863 0.8803149
## 22024 1.2292565 0.2036947 1.098644 0.010656575 0.3137328
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 22019 0.03301847 0.08246113 0.3444263 21
## 22020 0.02319654 0.10291531 0.4108362 21
## 22021 0.01345187 0.05699807 0.4067141 21
## 22022 0.04165140 0.10278206 0.3426506 21
## 22023 0.04035386 0.09667945 0.3381935 21
## 22024 0.03034663 0.10294394 0.3840779 21
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 31
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 33028 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 33029 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 33030 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 33031 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 33032 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 33033 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega
## 33028 2.0960974 0.3006922 1.198342 0.011623617 0.6564670
## 33029 2.7653413 0.3006922 1.269609 0.012314892 0.2751829
## 33030 1.1814994 0.3006922 1.107371 0.010741226 0.3211240
## 33031 0.2321777 0.3006922 1.020244 0.009896113 0.1560016
## 33032 1.5218109 0.3006922 1.140384 0.011061443 0.8875340
## 33033 1.2292565 0.3006922 1.111946 0.010785597 0.3290808
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 33028 0.03344512 0.08060069 0.3390899 31
## 33029 0.02353034 0.09936501 0.4064032 31
## 33030 0.01392642 0.05637268 0.4024378 31
## 33031 0.04161559 0.09924506 0.3374906 31
## 33032 0.04016319 0.09371839 0.3341247 31
## 33033 0.03057830 0.09939078 0.3791716 31
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 41
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 44037 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 44038 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 44039 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 44040 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 44041 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 44042 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega
## 44037 2.0960974 0.3976897 1.216126 0.011796117 0.6284494
## 44038 2.7653413 0.3976897 1.294525 0.012556573 0.2515237
## 44039 1.1814994 0.3976897 1.116605 0.010830788 0.2951232
## 44040 0.2321777 0.3976897 1.021910 0.009912274 0.1406015
## 44041 1.5218109 0.3976897 1.152647 0.011180385 0.8747665
## 44042 1.2292565 0.3976897 1.121594 0.010879181 0.3027231
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 44037 0.03269998 0.08382306 0.3482287 41
## 44038 0.02295863 0.10560023 0.4139441 41
## 44039 0.01311178 0.05744207 0.4097217 41
## 44040 0.04167657 0.10545630 0.3463468 41
## 44041 0.04050040 0.09889502 0.3410757 41
## 44042 0.03018045 0.10563117 0.3875481 41
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 51
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 55046 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 55047 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 55048 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 55049 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 55050 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 55051 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega
## 55046 2.0960974 0.4946872 1.230192 0.011932555 0.6332387
## 55047 2.7653413 0.4946872 1.314315 0.012748529 0.2554151
## 55048 1.1814994 0.4946872 1.123866 0.010901223 0.2994194
## 55049 0.2321777 0.4946872 1.023213 0.009924908 0.1431049
## 55050 1.5218109 0.4946872 1.162311 0.011274123 0.8770022
## 55051 1.2292565 0.4946872 1.129184 0.010952800 0.3070816
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 55046 0.03282735 0.08328107 0.3467258 51
## 55047 0.02305266 0.10452268 0.4127207 51
## 55048 0.01324638 0.05726675 0.4085369 51
## 55049 0.04166666 0.10438310 0.3448838 51
## 55050 0.04044135 0.09800831 0.3399382 51
## 55051 0.03024624 0.10455267 0.3861790 51by(new.result.df2, new.result.df2[, "rho.class"], summary)## new.result.df2[, "rho.class"]: 1
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.0097 Min. :1
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.0097 1st Qu.:1
## Median :0.21634 Median : 1.6936 Median :0.0097 Median :1
## Mean :0.29881 Mean : 2.9708 Mean :0.0097 Mean :1
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.0097 3rd Qu.:1
## Max. :1.00000 Max. :39.7812 Max. :0.0097 Max. :1
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.0097 Min. :0.09131 Min. :0.004582 Min. :0.02976
## 1st Qu.:0.0097 1st Qu.:0.36400 1st Qu.:0.024829 1st Qu.:0.07093
## Median :0.0097 Median :0.52523 Median :0.035162 Median :0.09819
## Mean :0.0097 Mean :0.54385 Mean :0.045850 Mean :0.08474
## 3rd Qu.:0.0097 3rd Qu.:0.70825 3rd Qu.:0.051346 3rd Qu.:0.09845
## Max. :0.0097 Max. :0.99758 Max. :1.223416 Max. :0.09849
##
## p.below.0 rho.class
## Min. :0.0000 1 :11009
## 1st Qu.:0.2707 11: 0
## Median :0.3346 21: 0
## Mean :0.3117 31: 0
## 3rd Qu.:0.3794 41: 0
## Max. :0.4779 51: 0
##
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 11
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.1067 Min. : 1.000
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.1067 1st Qu.: 1.050
## Median :0.21634 Median : 1.6936 Median :0.1067 Median : 1.107
## Mean :0.29881 Mean : 2.9708 Mean :0.1067 Mean : 1.233
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.1067 3rd Qu.: 1.238
## Max. :1.00000 Max. :39.7812 Max. :0.1067 Max. :11.000
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.00970 Min. :0.07839 Min. :0.004028 Min. :0.02997
## 1st Qu.:0.01018 1st Qu.:0.32635 1st Qu.:0.023889 1st Qu.:0.07394
## Median :0.01074 Median :0.48359 Median :0.034613 Median :0.10666
## Mean :0.01196 Mean :0.50969 Mean :0.045565 Mean :0.09073
## 3rd Qu.:0.01201 3rd Qu.:0.67266 3rd Qu.:0.051252 3rd Qu.:0.10700
## Max. :0.10670 Max. :0.99715 Max. :1.266254 Max. :0.10705
##
## p.below.0 rho.class
## Min. :0.0000 1 : 0
## 1st Qu.:0.2837 11:11009
## Median :0.3461 21: 0
## Mean :0.3232 31: 0
## 3rd Qu.:0.3905 41: 0
## Max. :0.4816 51: 0
##
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 21
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.2037 Min. : 1.000
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.2037 1st Qu.: 1.064
## Median :0.21634 Median : 1.6936 Median :0.2037 Median : 1.138
## Mean :0.29881 Mean : 2.9708 Mean :0.2037 Mean : 1.323
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.2037 3rd Qu.: 1.311
## Max. :1.00000 Max. :39.7812 Max. :0.2037 Max. :21.000
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.00970 Min. :0.08408 Min. :0.004274 Min. :0.02988
## 1st Qu.:0.01032 1st Qu.:0.34334 1st Qu.:0.024329 1st Qu.:0.07258
## Median :0.01104 Median :0.50265 Median :0.034861 Median :0.10270
## Mean :0.01283 Mean :0.52532 Mean :0.045688 Mean :0.08794
## 3rd Qu.:0.01272 3rd Qu.:0.68922 3rd Qu.:0.051270 3rd Qu.:0.10300
## Max. :0.20369 Max. :0.99736 Max. :1.246068 Max. :0.10304
##
## p.below.0 rho.class
## Min. :0.0000 1 : 0
## 1st Qu.:0.2778 11: 0
## Median :0.3409 21:11009
## Mean :0.3181 31: 0
## 3rd Qu.:0.3856 41: 0
## Max. :0.4800 51: 0
##
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 31
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.3007 Min. : 1.000
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.3007 1st Qu.: 1.072
## Median :0.21634 Median : 1.6936 Median :0.3007 Median : 1.157
## Mean :0.29881 Mean : 2.9708 Mean :0.3007 Mean : 1.386
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.3007 3rd Qu.: 1.358
## Max. :1.00000 Max. :39.7812 Max. :0.3007 Max. :31.000
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.00970 Min. :0.08966 Min. :0.004512 Min. :0.02979
## 1st Qu.:0.01040 1st Qu.:0.35938 1st Qu.:0.024718 1st Qu.:0.07130
## Median :0.01123 Median :0.52023 Median :0.035112 Median :0.09917
## Mean :0.01344 Mean :0.53974 Mean :0.045813 Mean :0.08544
## 3rd Qu.:0.01317 3rd Qu.:0.70410 3rd Qu.:0.051296 3rd Qu.:0.09944
## Max. :0.30069 Max. :0.99753 Max. :1.228319 Max. :0.09948
##
## p.below.0 rho.class
## Min. :0.0000 1 : 0
## 1st Qu.:0.2724 11: 0
## Median :0.3360 21: 0
## Mean :0.3132 31:11009
## 3rd Qu.:0.3809 41: 0
## Max. :0.4784 51: 0
##
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 41
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.3977 Min. : 1.000
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.3977 1st Qu.: 1.078
## Median :0.21634 Median : 1.6936 Median :0.3977 Median : 1.171
## Mean :0.29881 Mean : 2.9708 Mean :0.3977 Mean : 1.436
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.3977 3rd Qu.: 1.392
## Max. :1.00000 Max. :39.7812 Max. :0.3977 Max. :41.000
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.00970 Min. :0.08019 Min. :0.004106 Min. :0.02994
## 1st Qu.:0.01046 1st Qu.:0.33179 1st Qu.:0.024027 1st Qu.:0.07350
## Median :0.01136 Median :0.48974 Median :0.034684 Median :0.10537
## Mean :0.01393 Mean :0.51474 Mean :0.045603 Mean :0.08982
## 3rd Qu.:0.01350 3rd Qu.:0.67806 3rd Qu.:0.051232 3rd Qu.:0.10569
## Max. :0.39769 Max. :0.99722 Max. :1.259626 Max. :0.10574
##
## p.below.0 rho.class
## Min. :0.0000 1 : 0
## 1st Qu.:0.2818 11: 0
## Median :0.3445 21: 0
## Mean :0.3216 31: 0
## 3rd Qu.:0.3891 41:11009
## Max. :0.4811 51: 0
##
## --------------------------------------------------------
## new.result.df2[, "rho.class"]: 51
## ensembl_id betas.hat liver_residual_variance
## ENSMUSG00000000001: 1 Min. :0.00000 Min. :0.000975
## ENSMUSG00000000003: 1 1st Qu.:0.01186 1st Qu.:0.005551
## ENSMUSG00000000037: 1 Median :0.02430 Median :0.010731
## ENSMUSG00000000049: 1 Mean :0.05436 Mean :0.029009
## ENSMUSG00000000056: 1 3rd Qu.:0.05360 3rd Qu.:0.023547
## ENSMUSG00000000058: 1 Max. :2.76340 Max. :4.004622
## (Other) :11003
## liver_pvalue neg_log_lung_pvalue rho tmm
## Min. :0.00000 Min. : 0.0000 Min. :0.4947 Min. : 1.000
## 1st Qu.:0.03867 1st Qu.: 0.8062 1st Qu.:0.4947 1st Qu.: 1.083
## Median :0.21634 Median : 1.6936 Median :0.4947 Median : 1.182
## Mean :0.29881 Mean : 2.9708 Mean :0.4947 Mean : 1.478
## 3rd Qu.:0.50673 3rd Qu.: 3.5428 3rd Qu.:0.4947 3rd Qu.: 1.419
## Max. :1.00000 Max. :39.7812 Max. :0.4947 Max. :51.000
##
## tau omega beta_tieda n.betas.tieda.se
## Min. :0.00970 Min. :0.08172 Min. :0.004172 Min. :0.02992
## 1st Qu.:0.01050 1st Qu.:0.33637 1st Qu.:0.024136 1st Qu.:0.07314
## Median :0.01147 Median :0.49488 Median :0.034764 Median :0.10430
## Mean :0.01433 Mean :0.51895 Mean :0.045636 Mean :0.08907
## 3rd Qu.:0.01377 3rd Qu.:0.68253 3rd Qu.:0.051227 3rd Qu.:0.10461
## Max. :0.49469 Max. :0.99727 Max. :1.254172 Max. :0.10466
##
## p.below.0 rho.class
## Min. :0.0000 1 : 0
## 1st Qu.:0.2804 11: 0
## Median :0.3431 21: 0
## Mean :0.3202 31: 0
## 3rd Qu.:0.3878 41: 0
## Max. :0.4806 51:11009
## ```
qplot(omega, 1-p.below.0, group=rho.class, colour = rho.class, shape=rho.class, data = new.result.df)
qplot(omega, beta_tieda, group=rho.class, colour = rho.class, shape=rho.class, data = new.result.df)
qplot(betas.hat, beta_tieda, group=rho.class, colour = rho.class, shape=rho.class, data = new.result.df)
qplot(neg_log_lung_pvalue,tmm, group=rho.class, colour = rho.class, shape=rho.class, data = new.result.df)
qplot(neg_log_lung_pvalue, data = new.result.df)
qplot(exp(-neg_log_lung_pvalue), data = new.result.df)
qplot(tmm, group=rho.class, data = new.result.df)
sorted.new.result.df2 <- new.result.df2[order(new.result.df2$rho.class, new.result.df2$neg_log_lung_pvalue),]
by(sorted.new.result.df2[c(1, 2,5, 10)], sorted.new.result.df2[, "rho.class"], tail)## sorted.new.result.df2[, "rho.class"]: 1
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 4868 ENSMUSG00000028656 1.064626984 29.46768 0.7888386
## 8589 ENSMUSG00000044827 0.043061086 29.63010 0.1112840
## 8658 ENSMUSG00000045594 0.272539683 29.80771 0.1810885
## 1614 ENSMUSG00000020022 0.581488038 30.48351 0.2558550
## 993 ENSMUSG00000011884 0.003957014 32.50629 0.1357632
## 8387 ENSMUSG00000042684 0.058696429 39.78123 0.1519445
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 11
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 15877 ENSMUSG00000028656 1.064626984 29.46768 0.8050848
## 19598 ENSMUSG00000044827 0.043061086 29.63010 0.1030722
## 19667 ENSMUSG00000045594 0.272539683 29.80771 0.1849379
## 12623 ENSMUSG00000020022 0.581488038 30.48351 0.2645083
## 12002 ENSMUSG00000011884 0.003957014 32.50629 0.1287301
## 19396 ENSMUSG00000042684 0.058696429 39.78123 0.1418963
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 21
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 26886 ENSMUSG00000028656 1.064626984 29.46768 0.7975607
## 30607 ENSMUSG00000044827 0.043061086 29.63010 0.1067378
## 30676 ENSMUSG00000045594 0.272539683 29.80771 0.1831382
## 23632 ENSMUSG00000020022 0.581488038 30.48351 0.2604279
## 23011 ENSMUSG00000011884 0.003957014 32.50629 0.1319977
## 30405 ENSMUSG00000042684 0.058696429 39.78123 0.1464170
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 31
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 37895 ENSMUSG00000028656 1.064626984 29.46768 0.7907520
## 41616 ENSMUSG00000044827 0.043061086 29.63010 0.1102577
## 41685 ENSMUSG00000045594 0.272539683 29.80771 0.1815348
## 34641 ENSMUSG00000020022 0.581488038 30.48351 0.2568441
## 34020 ENSMUSG00000011884 0.003957014 32.50629 0.1349392
## 41414 ENSMUSG00000042684 0.058696429 39.78123 0.1507043
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 41
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 48904 ENSMUSG00000028656 1.064626984 29.46768 0.8026397
## 52625 ENSMUSG00000044827 0.043061086 29.63010 0.1042388
## 52694 ENSMUSG00000045594 0.272539683 29.80771 0.1843498
## 45650 ENSMUSG00000020022 0.581488038 30.48351 0.2631681
## 45029 ENSMUSG00000011884 0.003957014 32.50629 0.1297939
## 52423 ENSMUSG00000042684 0.058696429 39.78123 0.1433414
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 51
## ensembl_id betas.hat neg_log_lung_pvalue beta_tieda
## 59913 ENSMUSG00000028656 1.064626984 29.46768 0.8006089
## 63634 ENSMUSG00000044827 0.043061086 29.63010 0.1052256
## 63703 ENSMUSG00000045594 0.272539683 29.80771 0.1838637
## 56659 ENSMUSG00000020022 0.581488038 30.48351 0.2620655
## 56038 ENSMUSG00000011884 0.003957014 32.50629 0.1306761
## 63432 ENSMUSG00000042684 0.058696429 39.78123 0.1445590by(sorted.new.result.df2, sorted.new.result.df2[, "rho.class"], tail)## sorted.new.result.df2[, "rho.class"]: 1
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 4868 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 8589 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 8658 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 1614 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 993 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 8387 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega beta_tieda
## 4868 29.46768 0.009699749 1 0.009699749 0.6548994 0.7888386
## 8589 29.63010 0.009699749 1 0.009699749 0.2455859 0.1112840
## 8658 29.80771 0.009699749 1 0.009699749 0.7577397 0.1810885
## 1614 30.48351 0.009699749 1 0.009699749 0.8494937 0.2558550
## 993 32.50629 0.009699749 1 0.009699749 0.6892381 0.1357632
## 8387 39.78123 0.009699749 1 0.009699749 0.3341640 0.1519445
## n.betas.tieda.se p.below.0 rho.class
## 4868 0.09846466 5.671021e-16 1
## 8589 0.04880698 1.130136e-02 1
## 8658 0.09847356 3.296088e-02 1
## 1614 0.09847969 4.687769e-03 1
## 993 0.08176452 4.841523e-02 1
## 8387 0.05693248 3.805517e-03 1
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 11
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 15877 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 19598 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 19667 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 12623 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 12002 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 19396 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega
## 15877 29.46768 0.1066972 5.907512 0.05730139 0.6163203
## 19598 29.63010 0.1066972 5.965634 0.05786515 0.2160253
## 19667 29.80771 0.1066972 6.029842 0.05848796 0.7258445
## 12623 30.48351 0.1066972 6.280542 0.06091968 0.8269193
## 12002 32.50629 0.1066972 7.094953 0.06881926 0.6524606
## 19396 39.78123 0.1066972 11.000000 0.10669724 0.2981552
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 15877 0.8050848 0.10701824 2.679227e-14 11
## 19598 0.1030722 0.04975401 1.914972e-02 11
## 19667 0.1849379 0.10702967 4.200197e-02 11
## 12623 0.2645083 0.10703754 6.733524e-03 11
## 12002 0.1287301 0.08646751 6.827413e-02 11
## 19396 0.1418963 0.05845168 7.599903e-03 11
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 21
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 26886 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 30607 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 30676 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 23632 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 23011 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 30405 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega
## 26886 29.46768 0.2036947 9.537292 0.09250934 0.6341874
## 30607 29.63010 0.2036947 9.656586 0.09366647 0.2292207
## 30676 29.80771 0.2036947 9.788739 0.09494832 0.7407567
## 23632 30.48351 0.2036947 10.308336 0.09998828 0.8375640
## 23011 32.50629 0.2036947 12.034287 0.11672956 0.6695471
## 30405 39.78123 0.2036947 21.000000 0.20369474 0.3143557
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 26886 0.7975607 0.10301627 4.889969e-15 21
## 30607 0.1067378 0.04933352 1.524767e-02 21
## 30676 0.1831382 0.10302647 3.773614e-02 21
## 23632 0.2604279 0.10303349 5.742178e-03 21
## 23011 0.1319977 0.08431517 5.872972e-02 21
## 30405 0.1464170 0.05777312 5.632883e-03 21
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 31
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 37895 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 41616 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 41685 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 34641 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 34020 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 41414 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega
## 37895 29.46768 0.3006922 12.72671 0.1234459 0.6503558
## 41616 29.63010 0.3006922 12.90641 0.1251889 0.2418916
## 41685 29.80771 0.3006922 13.10580 0.1271230 0.7540416
## 34641 30.48351 0.3006922 13.89309 0.1347595 0.8469133
## 34020 32.50629 0.3006922 16.54364 0.1604692 0.6849291
## 41414 39.78123 0.3006922 31.00000 0.3006922 0.3297196
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 37895 0.7907520 0.09945587 9.266838e-16 31
## 41616 0.1102577 0.04892634 1.211238e-02 31
## 41685 0.1815348 0.09946504 3.399210e-02 31
## 34641 0.2568441 0.09947136 4.910183e-03 31
## 34020 0.1349392 0.08232943 5.060501e-02 31
## 41414 0.1507043 0.05712217 4.166380e-03 31
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 41
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 48904 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 52625 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 52694 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 45650 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 45029 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 52423 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega
## 48904 29.46768 0.3976897 15.65521 0.1518516 0.6221267
## 52625 29.63010 0.3976897 15.89439 0.1541716 0.2202250
## 52694 29.80771 0.3976897 16.16010 0.1567489 0.7307176
## 45650 30.48351 0.3976897 17.21242 0.1669562 0.8304155
## 45029 32.50629 0.3976897 20.78971 0.2016549 0.6580235
## 52423 39.78123 0.3976897 41.00000 0.3976897 0.3033338
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 48904 0.8026397 0.10570931 1.564178e-14 41
## 52625 0.1042388 0.04962057 1.783283e-02 41
## 52694 0.1843498 0.10572033 4.060140e-02 41
## 45650 0.2631681 0.10572791 6.403213e-03 41
## 45029 0.1297939 0.08577270 6.511044e-02 41
## 52423 0.1433414 0.05823563 6.919729e-03 41
## --------------------------------------------------------
## sorted.new.result.df2[, "rho.class"]: 51
## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 59913 ENSMUSG00000028656 1.064626984 0.018407266 4.552986e-26
## 63634 ENSMUSG00000044827 0.043061086 0.003157578 2.733728e-04
## 63703 ENSMUSG00000045594 0.272539683 0.030338790 1.481601e-08
## 56659 ENSMUSG00000020022 0.581488038 0.054747720 1.767689e-13
## 56038 ENSMUSG00000011884 0.003957014 0.021513047 8.846174e-01
## 63432 ENSMUSG00000042684 0.058696429 0.004868026 8.328277e-05
## neg_log_lung_pvalue rho tmm tau omega
## 59913 29.46768 0.4946872 18.40226 0.1784973 0.6269491
## 63634 29.63010 0.4946872 18.70006 0.1813859 0.2237770
## 63703 29.80771 0.4946872 19.03122 0.1845981 0.7347450
## 56659 30.48351 0.4946872 20.34580 0.1973492 0.8332920
## 56038 32.50629 0.4946872 24.84853 0.2410245 0.6626363
## 63432 39.78123 0.4946872 51.00000 0.4946872 0.3076973
## beta_tieda n.betas.tieda.se p.below.0 rho.class
## 59913 0.8006089 0.10462845 9.899718e-15 51
## 63634 0.1052256 0.04950742 1.677452e-02 51
## 63703 0.1838637 0.10463913 3.944853e-02 51
## 56659 0.2620655 0.10464648 6.134815e-03 51
## 56038 0.1306761 0.08519227 6.252764e-02 51
## 63432 0.1445590 0.05805297 6.384946e-03 51liver.ASE <- read.csv(file= "ASE.genetics.113.153882-6.csv")
liver.ASE.symbol <- unique(liver.ASE$geneID)
length(liver.ASE.symbol)## [1] 440#liver.ASE <-liver.ASE[liver.ASE$geneID %in% names(table(liver.ASE$geneID))[table(liver.ASE$geneID) >=2], ]
sub.liver.ASE <-liver.ASE[which(liver.ASE$pvalBH.DxB7 < 0.05 & liver.ASE$pvalBH.BxD7 < 0.05), ]
dim(sub.liver.ASE)## [1] 1044 19liver.ASE.symbol <- unique(liver.ASE$geneID)
liver.ASE.symbol <- unique(sub.liver.ASE$geneID)
length(liver.ASE.symbol)## [1] 367head(liver.ASE.symbol)## [1] Rrs1 6330578E17Rik Aox3 Cps1 Ugt1a10
## [6] Heatr7b1
## 440 Levels: 0610010D20Rik 1100001G20Rik 1100001H23Rik ... Zfp180library(biomaRt)
mouse = useMart("ensembl", dataset = "mmusculus_gene_ensembl")
liver.ASE.ensembl <- getBM( attributes=c("ensembl_gene_id", "mgi_symbol") , filters=
"mgi_symbol" , values =liver.ASE.symbol ,mart=mouse)
dim(liver.ASE.ensembl)## [1] 327 2length(unique(liver.ASE.ensembl$ensembl_gene_id))## [1] 327new.result.df2.rho1 <- new.result.df2[new.result.df2$rho.class == 1, ]
head(new.result.df2.rho1)## ensembl_id betas.hat liver_residual_variance liver_pvalue
## 1 ENSMUSG00000000001 0.015986111 0.018910763 5.377718e-01
## 2 ENSMUSG00000000003 0.016880682 0.003757142 1.929225e-01
## 3 ENSMUSG00000000037 0.003865079 0.004681090 7.788138e-01
## 4 ENSMUSG00000000049 0.042233333 0.001829162 9.084753e-06
## 5 ENSMUSG00000000056 0.063604072 0.078095912 2.269175e-01
## 6 ENSMUSG00000000058 0.025611111 0.004853968 5.846922e-02
## neg_log_lung_pvalue rho tmm tau omega beta_tieda
## 1 2.0960974 0.009699749 1 0.009699749 0.6609725 0.03356494
## 2 2.7653413 0.009699749 1 0.009699749 0.2791983 0.02362737
## 3 1.1814994 0.009699749 1 0.009699749 0.3255088 0.01406380
## 4 0.2321777 0.009699749 1 0.009699749 0.1586587 0.04160507
## 5 1.5218109 0.009699749 1 0.009699749 0.8895190 0.04011077
## 6 1.2292565 0.009699749 1 0.009699749 0.3335208 0.03064532
## n.betas.tieda.se p.below.0 rho.class
## 1 0.08007039 0.3375376 1
## 2 0.09837648 0.4050981 1
## 3 0.05619033 0.4011821 1
## 4 0.09826008 0.3359952 1
## 5 0.09288763 0.3329360 1
## 6 0.09840149 0.3777363 1new.result.df2.rho1$rank.p.below.0 <- rank(new.result.df2.rho1$p.below.0)
new.result.df2.rho1$rank.liver.pvalue <- rank(new.result.df2.rho1$liver_pvalue)
result.rho1 <- matrix(, nrow(new.result.df2.rho1), 8)
colnames(result.rho1)<-c("bayrank","bayover","bay_TPR","bay_FPR", "orirank","oriover","ori_TPR","ori_FPR" )
for (i in 1:nrow(new.result.df2.rho1))
{
newdata1.rho1 <- subset(new.result.df2.rho1, rank.p.below.0 <= i)
overlap.newdata1.rho1 <- newdata1.rho1[newdata1.rho1$ensembl_id %in% liver.ASE.ensembl$ensembl_gene_id, ]
result.rho1[i, 1] <- i
result.rho1[i, 2] <- nrow(overlap.newdata1.rho1)/nrow(newdata1.rho1)
result.rho1[i, 3] <- nrow(overlap.newdata1.rho1)/nrow(liver.ASE.ensembl)
newdata2.rho1 <- subset(new.result.df2.rho1, rank.liver.pvalue <= i)
overlap.newdata2.rho1 <- newdata2.rho1[newdata2.rho1$ensembl_id %in% liver.ASE.ensembl$ensembl_gene_id, ]
result.rho1[i, 5] <- i
result.rho1[i, 6] <- nrow(overlap.newdata2.rho1)/nrow(newdata2.rho1)
result.rho1[i, 7] <- nrow(overlap.newdata2.rho1)/nrow(liver.ASE.ensembl)
}
head(result.rho1)## bayrank bayover bay_TPR bay_FPR orirank oriover ori_TPR
## [1,] 1 0.0000000 0.000000000 NA 1 0.0000000 0.000000000
## [2,] 2 0.5000000 0.003058104 NA 2 0.5000000 0.003058104
## [3,] 3 0.3333333 0.003058104 NA 3 0.6666667 0.006116208
## [4,] 4 0.2500000 0.003058104 NA 4 0.5000000 0.006116208
## [5,] 5 0.4000000 0.006116208 NA 5 0.4000000 0.006116208
## [6,] 6 0.3333333 0.006116208 NA 6 0.5000000 0.009174312
## ori_FPR
## [1,] NA
## [2,] NA
## [3,] NA
## [4,] NA
## [5,] NA
## [6,] NAtail(result.rho1)## bayrank bayover bay_TPR bay_FPR orirank oriover ori_TPR
## [11004,] 11004 0.02417303 0.8134557 NA 11004 0.02417303 0.8134557
## [11005,] 11005 0.02417083 0.8134557 NA 11005 0.02417083 0.8134557
## [11006,] 11006 0.02416864 0.8134557 NA 11006 0.02416864 0.8134557
## [11007,] 11007 0.02416644 0.8134557 NA 11007 0.02416644 0.8134557
## [11008,] 11008 0.02416424 0.8134557 NA 11008 0.02416424 0.8134557
## [11009,] 11009 0.02416205 0.8134557 NA 11009 0.02416205 0.8134557
## ori_FPR
## [11004,] NA
## [11005,] NA
## [11006,] NA
## [11007,] NA
## [11008,] NA
## [11009,] NAplot(result.rho1[, 1], result.rho1[, 3], type="l", col="red", xlab="Ranking", ylab="TPR", ylim=c(0, 1), xlim=c(0, 1800) )
par(new=TRUE)
plot( result.rho1[, 1], result.rho1[, 7], type="l", col="green", xlab="Ranking", ylab="TPR", ylim=c(0, 1), xlim=c(0, 1800) )
legend("topright", inset=.05, c("Bayesian","Original"), text.col = c("red", "green"), horiz=TRUE)
plot(result[, 1], result[, 2], type=“l”, col=“red”, xlab=“Ranking”, ylab=“Percent of true positive”, ylim=c(0, 1))
par(new=TRUE)
plot(result[, 5], result[, 6], type=“l”, col=“green”, xlab=“Ranking”, ylab=“Percent of true positive”, ylim=c(0, 1))
legend(“topright”, inset=.05, c(“Bayesian”,“Original”), text.col = c(“red”, “green”), horiz=TRUE)