balliFit function

balliFit

Estimates likelihood and Bartlett correction factor using BALLI algorithm of each gene

balliFit(y_mat, x_mat, tecVar, intVar = 2, full = T, cfault = 0,
  miter = 200, conv = 1e-06)

Arguments

• y_mat: numeric vector containing log-cpm values of each gene and each sample
• x_mat: design matrix with samples in row and covariable(s) to be estimated in column
• tecVar: numeric vector containing estimated technical variance of a gene of each sample
• intVar: numeric vector designating interest variable(s) which is(are) column number(s) of x_mat
• full: logical value designating full model (TRUE) or reduced model (FALSE).
• cfault: initial value of index showing whether converged (0) or not (1).
• miter: maximum number of iteration to converge.
• conv: threshold for convergence

Returns

following components are estimated - ll: log-likelihoods

• beta: coefficients of interested variable(s)

• alpha: coefficients of nuisance variable(s)

• BCF: Bartlett's correction factor

• cfault: index whether converged or not

Examples

expr <- data.frame(t(sapply(1:1000,function(x)rnbinom(20,mu=500,size=50))))
group <- c(rep("A",10),rep("B",10))
design <- model.matrix(~group, data = expr)
dge <- DGEList(counts=expr, group=group)
dge <- calcNormFactors(dge)
tV <- tecVarEstim(dge,design)
gtv <- tV$tecVar[1,]
gdat <- data.frame(logcpm=tV$logcpm[1,],design,tecVar=gtv)
gy <- matrix(unlist(gdat[,1]),ncol=1)
gx <- matrix(unlist(gdat[,2:(ncol(gdat)-1)]),ncol=ncol(gdat)-2)
balliFit(y_mat=gy,x_mat=gx,tecVar=gtv,intVar=2,full=TRUE,cfault=0,miter=200,conv=1e-6)