QuantileFirstPrinciples.R

From Organic Design wiki
Revision as of 07:27, 29 May 2007 by Nad (talk | contribs) ({{R}})
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
  1. ===============================================================================
  2. Generating Example Data for quantile normalization

Code snipits and programs written in R, S or S-PLUS

  1. ===============================================================================

y <- cbind(ProbeSet = rep(1:2, each=3), Probe = rep(1:3, 2),

          Chip1 = c(7,3,2,4,10,12), Chip2 = c(9,5,6,8,11,10),
          Chip3 = c(19,14,11,8,16,15))

library(limma) ynorm <- cbind(y[,1:2],normalizeQuantiles(y[,3:5], ties=T))

  1. ===============================================================================
  2. Generating Example Data for median polish
  3. ===============================================================================

y <- cbind(probe1=c(18,13,15,19), probe2=c(11,7,6,15), probe3=c(8,5,7,12),

          probe4=c(21,16,16,18), probe5=c(4,7,6,5))

E <- y

  1. ===============================================================================
  2. Median polish by first principles (3x), convergence is fast
  3. ===============================================================================

FUN <- median

  1. Sweep out medians on rows

meds <- apply(E, 1, FUN) E <- sweep(E, 1, STATS=meds)

  1. Sweep out medians on columns

meds <- apply(E, 2, FUN) E <- sweep(E, 2, STATS=meds)

  1. Sweep out medians on rows

meds <- apply(E, 1, FUN) E <- sweep(E, 1, STATS=meds)

  1. Sweep out medians on columns

meds <- apply(E, 2, FUN) E <- sweep(E, 2, STATS=meds)

  1. Sweep out medians on rows

meds <- apply(E, 1, FUN) E <- sweep(E, 1, STATS=meds)

  1. Sweep out medians on columns

meds <- apply(E, 2, FUN) E <- sweep(E, 2, STATS=meds)

  1. see ?medpolish for R function

medpolish(y, trace.iter=F)$res

  1. ===============================================================================
  2. Fitted values
  3. ===============================================================================

mu <- y - E

  1. ===============================================================================
  2. Y_ijk = mu_i.k + alpha_.jk + Eijk -Gene chip specific probe affinities mu_i.k
  3. i = 1..I chips, j = 1..J probes, k = 1..G genes
  4. ===============================================================================
  1. mu_i.k calculation

mu.i <- apply(mu,1, mean)

  1. alpha_.jk calculation

alpha <- y - mu.i - E alpha.j <- apply(alpha,2, mean)


  1. ===============================================================================
  2. Mean instead of median decomposition using anova, aov
  3. ===============================================================================

dset <- as.data.frame(cbind(y = c(y), probe = rep(1:5, each=4), array = rep(1:4, 5))) dset[,"probe"] <- factor(dset[,"probe"]) dset[,"array"] <- factor(dset[,"array"])

fit <- aov(y ~ probe + array, data=dset) summary(fit)

matrix(fit$resid, nc=5)