Difference between revisions of "LNNB-snippet.R"
From Organic Design wiki
m |
|||
| (2 intermediate revisions by the same user not shown) | |||
| Line 1: | Line 1: | ||
| − | p | + | # {{R}} |
| − | n | + | # p (genes) by n (slides) matrix |
| − | + | # m=p keeping fdr notation | |
| − | |||
| − | |||
| − | + | m <- p <- 1000 | |
| − | + | nreps <- 5 | |
| − | X | + | |
| − | dim(X) <- c(p | + | n <- nreps * 2 |
| − | X | + | pi0 <- 0.95 |
| + | |||
| + | m1 <- round(m * (1-pi0)) | ||
| + | m0 <- p - m1 | ||
| + | |||
| + | params <- (2.32501, sqrt(0.38106), sqrt(1.11561)) # "IPTG-a parameters" | ||
| + | k <- 1000 | ||
| + | |||
| + | mu <- params[1] | ||
| + | sigma <- params[2] | ||
| + | tau <- params[3] | ||
| + | |||
| + | DEmeans <- rnorm(2*m1, mu, tau) | ||
| + | EEmeans <- rnorm(m0 , mu, tau) | ||
| + | means <- c(rep(DEmeans, each=nreps),rep(EEmeans, each=2*nreps)) | ||
| + | X <- exp(rnorm(n*p, means, sigma)) | ||
| + | dim(X) <- c(n,p) | ||
| + | X <- t(X) | ||
| + | |||
| + | # Graphical check | ||
| + | hist(kDE, prob=TRUE) | ||
| + | lines(density(kDE)) | ||
| + | |||
| + | hist(kEE, prob=TRUE) | ||
| + | lines(density(kEE)) | ||
| + | |||
| + | Xbar <- rowMeans(log2(X)) | ||
| + | hist(Xbar, breaks=30, prob=TRUE) | ||
| + | lines(density(Xbar)) | ||
| + | pairs(log2(X)) | ||
Latest revision as of 00:39, 6 June 2007
Code snipits and programs written in R, S or S-PLUS
- p (genes) by n (slides) matrix
- m=p keeping fdr notation
m <- p <- 1000 nreps <- 5
n <- nreps * 2 pi0 <- 0.95
m1 <- round(m * (1-pi0)) m0 <- p - m1
params <- (2.32501, sqrt(0.38106), sqrt(1.11561)) # "IPTG-a parameters" k <- 1000
mu <- params[1] sigma <- params[2] tau <- params[3]
DEmeans <- rnorm(2*m1, mu, tau) EEmeans <- rnorm(m0 , mu, tau) means <- c(rep(DEmeans, each=nreps),rep(EEmeans, each=2*nreps)) X <- exp(rnorm(n*p, means, sigma)) dim(X) <- c(n,p) X <- t(X)
- Graphical check
hist(kDE, prob=TRUE) lines(density(kDE))
hist(kEE, prob=TRUE) lines(density(kEE))
Xbar <- rowMeans(log2(X)) hist(Xbar, breaks=30, prob=TRUE) lines(density(Xbar)) pairs(log2(X))
