### R code from vignette source 'vignettes/lpNet/inst/doc/vignette_lpNet.Rnw' ################################################### ### code chunk number 1: setup1 ################################################### library("lpNet") ################################################### ### code chunk number 2: data1 ################################################### n <- 5 # number of genes K <- 7 # number of perturbation experiments obs <- matrix(rnorm(n*K),nrow=n,ncol=K) delta <- apply(obs,1,mean) # perturbation vector (0 if gene inactivated and 1 otherwise) b <- c(0,1,1,1,1, # perturbation exp1 1,0,1,1,1, # perturbation exp2 1,1,0,1,1, # perturbation exp3... 1,1,1,0,1, 1,1,1,1,0, 1,0,0,1,1, 1,1,1,1,1) ################################################### ### code chunk number 3: lp ################################################### res1 <- doILP(obs,delta,lambda=1,b,n,K,annot=getEdgeAnnot(n)) adja1 <- getAdja(res1,n) ################################################### ### code chunk number 4: lp_pos ################################################### res2 <- doILP(obs,delta,lambda=1,b,n,K, annot=getEdgeAnnot(n,allpos=TRUE),all.pos=TRUE) adja2 <- getAdja(res2,n) ################################################### ### code chunk number 5: gaussian ################################################### active_mu <- 1.1 inactive_mu <- 0.9 active_sd <- inactive_sd <- 0.01 ################################################### ### code chunk number 6: loocv ################################################### times <- 5 # number of times the removed observation is sampled annot_node <- c(LETTERS[1:n]) # annotation of the nodes loocv_res <- loocv(times,obs,n,b,K,delta,lambda=1, annot=getEdgeAnnot(n),annot_node, active_mu,active_sd,inactive_mu,inactive_sd) loocv_res$MSE ################################################### ### code chunk number 7: loocv ################################################### adja_loocv <- getSampleAdjaMAD(loocv_res$edges_all,n,annot_node) ################################################### ### code chunk number 8: rangeLambda ################################################### lambda <- calcRangeLambda(delta,obs,stepsize=2) MSE <- Inf for(lamd in lambda){ loocv_res <- loocv(times,obs,n,b,K,delta,lambda=lamd, annot=getEdgeAnnot(n),annot_node, active_mu,active_sd,inactive_mu,inactive_sd) if(loocv_res$MSE",1)) res4 <- doILP(obs,delta,lambda=1,b,n,K, annot=getEdgeAnnot(n),prior=prior) adja4 <- getAdja(res4,n) ################################################### ### code chunk number 12: prior3 ################################################### prior <- list(c("w+_1_2",1/0.9,">",1)) res5 <- doILP(obs,delta,lambda=1,b,n,K, annot=getEdgeAnnot(n),prior=prior) adja5 <- getAdja(res5,n) ################################################### ### code chunk number 13: sahinData ################################################### data("SahinRNAi2008") dataStim <- dat.normalized[dat.normalized[,17]==1,-17] dataUnstim <- dat.normalized[dat.normalized[,17]==0,-17] # summarize replicates; transpose: rows=genes, cols=experiments dataSt <- t(summarizeRepl(dataStim,type=mean)) dataUnst <- t(summarizeRepl(dataUnstim,type=mean)) ################################################### ### code chunk number 14: sahinData_parameter ################################################### n <- 16 # number of genes K <- 16 # number of experiments annot <- getEdgeAnnot(n) # perturbation vector; kd of: b <- c(0,rep(1,15), # erbb1 0,0,rep(1,14), # erbb1 & 2 0,rep(1,14),0, # erbb1 & 3 1,0,rep(1,13),0, # erbb2 & 3 rep(1,10),0,1,1,1,1,1, # IGFR rep(1,11),0,1,1,1,1, # ERalpha rep(1,12),0,1,1,1, # MYC rep(1,7),0,rep(1,8), # AKT1 rep(1,8),0,rep(1,7), # MEK1 rep(1,5),0,rep(1,10), # CDK2 rep(1,6),0,rep(1,9), # CDK4 rep(1,13),0,1,1, # CDK6 1,1,0,rep(1,13), # p21 1,1,1,0,rep(1,12), # p27 rep(1,4),0,rep(1,11), # Cyclin D1 rep(1,14),0,1) # Cyclin E1 ################################################### ### code chunk number 15: sahinData_delta ################################################### delta <- as.double(dataUnst[,1]) delta[11:16] <- mean(dataUnst[,1],na.rm=T) ################################################### ### code chunk number 16: sahinData_lp ################################################### resERBB <- doILP(dataSt[,-1],delta,lambda=1.83,b,n,K,annot, all.pos=FALSE,sourceNode=c(1,2,16),sinkNode=10) adjaERBB <- getAdja(resERBB,n)