Nothing
.pls.loads <-
function(X, Y.lvs, blocks)
{
lvs <- length(blocks)
mvs <- ncol(X)
blocklist <- as.list(1:lvs)
for (j in 1:lvs)
blocklist[[j]] <- rep(j,blocks[j])
blocklist <- unlist(blocklist)
loads <- rep(NA, mvs)
comu <- rep(NA, mvs)
for (j in 1:lvs)
loads[blocklist==j] <- cor(X[,blocklist==j], Y.lvs[,j])
comu <- loads^2
names(loads) <- colnames(X)
names(comu) <- colnames(X)
res.loads <- list(loads, comu)
return(res.loads)
}
.pls.paths <-
function(IDM, Y.lvs, plsr)
{
lvs.names <- colnames(IDM)
endo = rowSums(IDM)
endo[endo!=0] <- 1 # vector indicating endogenous LVs
innmod <- as.list(1:sum(endo))
Path <- IDM
residuals <- as.list(1:sum(endo))
R2 <- rep(0,nrow(IDM))
for (aux in 1:sum(endo))
{
k1 <- which(endo==1)[aux] # index for endo LV
k2 <- which(IDM[k1,]==1) # index for indep LVs
if (length(k2)>1 & plsr) {
path.lm <- .plsr1(Y.lvs[,k2], Y.lvs[,k1], nc=2)
Path[k1,k2] <- path.lm$coeffs
residuals[[aux]] <- path.lm$resid
R2[k1] <- path.lm$R2[1]
inn.val <- c(path.lm$R2[1], path.lm$cte, path.lm$coeffs)
inn.lab <- c("R2", "Intercept", paste(rep("path_",length(k2)),names(k2),sep=""))
names(inn.val) <- NULL
innmod[[aux]] <- data.frame(concept=inn.lab, value=round(inn.val,4))
}
if (length(k2)==1 | !plsr) {
path.lm <- summary(lm(Y.lvs[,k1] ~ Y.lvs[,k2]))
Path[k1,k2] <- path.lm$coef[-1,1]
residuals[[aux]] <- path.lm$residuals
R2[k1] <- path.lm$r.squared
inn.val <- c(path.lm$r.squared, path.lm$coef[,1])
inn.lab <- c("R2", "Intercept", paste(rep("path_",length(k2)),names(k2),sep=""))
names(inn.val) <- NULL
innmod[[aux]] <- data.frame(concept=inn.lab, value=round(inn.val,4))
}
}
names(innmod) <- lvs.names[endo!=0]
names(R2) <- lvs.names
res.paths <- list(innmod, Path, R2, residuals)
return(res.paths)
}
.pls.weights <-
function(X, IDM, blocks, modes, scheme, tol, iter)
{
lvs <- nrow(IDM)
mvs <- ncol(X)
sdv <- sqrt((nrow(X)-1)/nrow(X)) # std.dev factor correction
blocklist <- as.list(1:lvs)
for (j in 1:lvs)
blocklist[[j]] <- rep(j,blocks[j])
blocklist <- unlist(blocklist)
# outer design matrix 'ODM' and matrix of outer weights 'W'
ODM <- matrix(0, mvs, lvs)
for (j in 1:lvs)
ODM[which(blocklist==j),j] <- rep(1,blocks[j])
W <- ODM %*% diag(1/(sd(X %*% ODM)*sdv),lvs,lvs)
w.old <- rowSums(W)
w.dif <- 1
itermax <- 1
repeat
{
Y <- X %*% W # external estimation of LVs 'Y'
Y <- scale(Y) * sdv
# matrix of inner weights 'e'
E <- switch(scheme,
"centroid" = sign(cor(Y) * (IDM + t(IDM))),
"factor" = cor(Y) * (IDM + t(IDM)),
"path" = .path.scheme(IDM, Y))
Z <- Y %*% E # internal estimation of LVs 'Z'
# scaling Z
Z <- Z %*% diag(1/(sd(Z)*sdv), lvs, lvs)
# computing outer weights 'w'
for (j in 1:lvs)
{
X.blok = X[,which(blocklist==j)]
if (modes[j]=="A")# reflective way
ODM[which(blocklist==j),j] <- (1/nrow(X)) * Z[,j] %*% X.blok
if (modes[j]=="B")# formative way
ODM[which(blocklist==j),j] <- solve.qr(qr(X.blok),Z[,j])
}
W <- ODM
w.new <- rowSums(W)
w.dif <- sum((abs(w.old) - abs(w.new))^2) # difference of out.weights
if (w.dif<tol || itermax==iter) break
w.old <- w.new
itermax <- itermax + 1
} # end repeat
W <- ODM %*% diag(1/(sd(X %*% ODM)*sdv),lvs,lvs)
w.new <- rowSums(W)
names(w.new) <- colnames(X)
dimnames(W) <- list(colnames(X),rownames(IDM))
res.ws <- list(w.new, W, itermax)
if (itermax==iter) res.ws=NULL
return(res.ws)
}
.plsr1 <-
function(x, y, nc=NULL, scaled=TRUE)
{
# ============ checking arguments ============
X <- as.matrix(x)
Y <- as.matrix(y)
n <- nrow(X)
p <- ncol(X)
if (is.null(nc))
nc <- p
# ============ setting inputs ==============
if (scaled) Xx<-scale(X) else Xx<-scale(X,scale=F)
if (scaled) Yy<-scale(Y) else Yy<-scale(Y,scale=F)
X.old <- Xx
Y.old <- Yy
Th <- matrix(NA, n, nc)# matrix of X-scores
Ph <- matrix(NA, p, nc)# matrix of X-loadings
Wh <- matrix(NA, p, nc)# matrix of raw-weights
Uh <- matrix(NA, n, nc)# matrix of Y-scores
ch <- rep(NA, nc)# vector of y-loadings
# ============ pls regression algorithm ==============
for (h in 1:nc)
{
w.old <- t(X.old) %*% Y.old / sum(Y.old^2)
w.new <- w.old / sqrt(sum(w.old^2)) # normalization
t.new <- X.old %*% w.new
p.new <- t(X.old) %*% t.new / sum(t.new^2)
c.new <- t(Y.old) %*% t.new / sum(t.new^2)
u.new <- Y.old / as.vector(c.new)
Y.old <- Y.old - t.new%*%c.new# deflate y.old
X.old <- X.old - (t.new %*% t(p.new))# deflate X.old
Th[,h] <- t.new
Ph[,h] <- p.new
Wh[,h] <- w.new
Uh[,h] <- u.new
ch[h] <- c.new
}
Ws <- Wh %*% solve(t(Ph)%*%Wh)# modified weights
Bs <- as.vector(Ws %*% ch) # std beta coeffs
Br <- Bs * (rep(sd(Y),p)/apply(X,2,sd)) # beta coeffs
cte <- as.vector(mean(y) - Br%*%apply(X,2,mean))# intercept
y.hat <- X%*%Br+cte# y predicted
resid <- as.vector(Y - y.hat)# residuals
R2 <- as.vector(cor(Th, Yy))^2 # R2 coefficients
names(Br) <- colnames(X)
names(resid) <- rownames(Y)
names(y.hat) <- rownames(Y)
names(R2) <- paste(rep("t",nc),1:nc,sep="")
res <- list(coeffs=Br, coef.std=Bs, cte=cte, R2=R2[1:nc], resid=resid, y.pred=y.hat)
return(res)
}
.path.scheme <-
function(IDM, Y)
{
lvs <- nrow(IDM)
E <- IDM
for (k in 1:lvs)
{
if (length(which(IDM[k,]==1)) > 0)
E[which(IDM[k,]==1),k] <- lm(Y[,k]~Y[,which(IDM[k,]==1)]-1)$coef
if (length(which(IDM[,k]==1)) > 0)
E[which(IDM[,k]==1),k] <- cor(Y[,k], Y[,which(IDM[,k]==1)])
}
return(E)
}
.nominal.split <-
function(v)
{
# INPUT(S)
# v: vector containing the nominal categories
# PARAM(S)
k = length(v) # length of vector
num.parts = (2^(k-1))-1 # number of binary partitions
# stop limit of order of partitions
if ((k%%2)==0) stop.parts=k/2 else stop.parts=(k-1)/2
part = as.list(1:num.parts) # list for warehousing the partitions
antipar = part # complementary list for part
# FUNCTION(S)
# internal function for evaluating candidates splits
eval.cand = function(x, y, aux)
{
# x: list of partitions
# y: candidate split to be included in part
# aux: num of parts already included in part
acum = 0
for (i in 1:(aux-1))
{
true.false = setequal(x[[i]], y)
if (true.false) acum=acum+1 else acum=acum+0
}
return(acum)
}
# first order partitions
aux = 1
for (i in 1:k)
{
part[[aux]] = v[i]
antipar[[aux]] = setdiff(v, part[[i]])
aux = aux + 1
}
# partitions of order greater than one
for (ord in 2:stop.parts)
{
if (stop.parts < ord)
break
# partitions of the immediate past order
p1 = aux - choose(k, (ord-1))
p2 = aux - 1
part.ant = part[p1:p2]
antipar.ant = antipar[p1:p2]
for (i in 1:length(part.ant))
{
part1 = part.ant[[i]]
anti1 = antipar.ant[[i]]
for (j in 1:length(anti1))
{
candidat = c(part1, anti1[j])
cand = eval.cand(part, candidat, aux)
if (cand==0)
{
part[[aux]] = candidat
antipar[[aux]] = setdiff(v, candidat)
aux = aux + 1
} else
{
next
}
}
}
}
# in case that k is even, only select partitions until num.parts
if (length(part) > num.parts)
{
part = part[1:num.parts]
antipar = antipar[1:num.parts]
}
parts = list(par1=part, par2=antipar)
return(parts)
}
.ordinal.split <-
function(cat.exe)
{
a = length(cat.exe)
if (a == 2) # binary variable
{
part = cat.exe[1]
antipar = cat.exe[2]
}
if (a > 2) # more than 2 categories
{
part = as.list(1:(a-1))
antipar = part
for (i in 1:(a-1))
{
part[[i]] = cat.exe[1:i]
antipar[[i]] = setdiff(cat.exe, part[[i]])
}
}
parts = list(par1=part, par2=antipar)
return(parts)
}
.pls.basic <-
function(DT, IDM, blocks, modes, scheme, scaled, tol, iter)
{
##### PARAMETERS
DM <- DT
lvs <- nrow(IDM)
mvs <- sum(blocks)
blocklist <- as.list(1:lvs)
for (j in 1:lvs)
blocklist[[j]] <- rep(j,blocks[j])
blocklist <- unlist(blocklist)
endo = rowSums(IDM)
endo[endo!=0] = 1
plsr <- FALSE
### Variable Names
lvs.names = colnames(IDM)
mvs.names = colnames(DM)
# apply the selected scaling
if (scaled) {
sd.X <- sqrt((nrow(DM)-1)/nrow(DM)) * apply(DM, 2, sd)
X <- scale(DM, scale=sd.X)
} else {
X <- scale(DM, scale=FALSE)
}
dimnames(X) <- list(rownames(DM), mvs.names)
# ==================== Stage 1: Iterative procedure ==================
out.ws <- .pls.weights(X, IDM, blocks, modes, scheme, tol, iter)
if (is.null(out.ws)) stop("The pls algorithm is non convergent")
out.weights <- round(out.ws[[1]], 4)
cor.XY <- cor(X, X%*%out.ws[[2]])
w.sig <- rep(NA,lvs)
for (k in 1:lvs)
w.sig[k] <- ifelse(sum(sign(cor.XY[which(blocklist==k),k]))<=0,-1,1)
Z.lvs <- X %*% out.ws[[2]] %*% diag(w.sig,lvs,lvs)
Y.lvs <- Z.lvs
if (!scaled)
Y.lvs <- DM %*% out.ws[[2]] %*% diag(w.sig,lvs,lvs)
dimnames(Y.lvs) <- list(rownames(X), lvs.names)
dimnames(Z.lvs) <- list(rownames(X), lvs.names)
loadcomu <- .pls.loads(X, Y.lvs, blocks)
loads <- loadcomu[[1]]
# ============ Stage 2: Path coefficients and total effects ==========
pathmod <- .pls.paths(IDM, Y.lvs, plsr)
innmod <- pathmod[[1]]
Path <- pathmod[[2]]
R2 <- pathmod[[3]]
residuals <- pathmod[[4]]
# ============================= Results ==============================
model <- list(IDM=IDM, blocks=blocks, scheme=scheme, modes=modes,
scaled=scaled, tol=tol, iter=iter)
resul = list(out.weights=out.weights, loadings=loads, scores=Y.lvs,
path.coefs=Path, R2=R2, residuals=residuals, model=model)
resul
}
.sanchez.aluja <-
function(Y.lvs, IDM, split.a, split.b)
{
lvs <- ncol(Y.lvs)
endo <- rowSums(IDM)
endo[endo!=0] <- 1
n.endo <- sum(endo)
n <- nrow(Y.lvs)
p.values <- rep(0,n.endo)
f.values <- rep(0,n.endo)
df1.values <- rep(0,n.endo)
df2.values <- rep(0,n.endo)
s1s2.values <- rep(0,n.endo)
who.endo <- which(endo==1)
SCE0 <- 0
SCE1 <- 0
Sa <- 0
Sb <- 0
p.sum <- 0
for (i in 1:n.endo)
{
aux <- who.endo[i]
indep <- which(IDM[aux,1:aux]==1)
reg <- cbind(Y.lvs[,aux], Y.lvs[,indep])
# SS0 calculation
lm.res <- lm(Y.lvs[,aux] ~ Y.lvs[,indep])$residuals
SS0 <- sum(lm.res^2)
SCE0 <- SCE0 + SS0
# regression matrices "Xa" and "Xb"
Xa <- reg[split.a,]
Xb <- reg[split.b,]
lm.resa <- lm(Xa[,1] ~ Xa[,2:ncol(Xa)])$residuals
lm.resb <- lm(Xb[,1] ~ Xb[,2:ncol(Xb)])$residuals
SS1 <- sum(lm.resa^2) + sum(lm.resb^2)
SCE1 <- SCE1 + SS1
Sa <- Sa + lm.resa
Sb <- Sb + lm.resb
# test
p <- length(indep)
p.sum <- p.sum + p
df1 <- p
df2 <- n - 2*p
f <- ((n-2*p)/p) * ((SS0-SS1)/SS1)
p.values[i] <- 1 - pf(f, df1, df2)
f.values[i] <- f
df1.values[i] <- df1
df2.values[i] <- df2
s1s2.values[i] <- round(sd(lm.resa)/sd(lm.resb),4)
}
n.sum <- n.endo * (length(split.a) + length(split.b))
f.sum <- ((n.sum-2*p.sum)/p.sum) * ((SCE0 - SCE1)/SCE1) # f global
df1.sum <- p.sum # gl num
df2.sum <- n.sum - 2*p.sum # gl den
pval.sum <- 1 - pf(f.sum, df1.sum, df2.sum) # pv
s1s2.sum <- sd(Sa) / sd(Sb)
F.global <- c(f.sum, df1.sum, df2.sum, pval.sum, s1s2.sum)
F.partial <- cbind(f.values, df1.values, df2.values, p.values, s1s2.values)
rownames(F.partial) <- rownames(IDM)[who.endo]
colnames(F.partial) <- c("f.stat", "df.num", "df.den", "p.val", "sa/sb")
res <- list(F.global, F.partial)
return(res)
}
.xexeloa <-
function(pls, DT, EXEV, type.exev, elemnod, nv, size, mox)
{
# cuauhmaitl elems
N <- nrow(EXEV)
if (nv==0) elems<-1:length(elemnod) else elems<-which(elemnod==nv)
indelem <- cbind(1:length(elems), elems)# identificador
E.node <- EXEV[elems,]# elems en cuauhmaitl actual
exevs <- apply(E.node, 2, function(x) nlevels(as.factor(x)))
list.cata <- as.list(1:length(exevs)) # categs a.cuauhmaitl
list.catb <- as.list(1:length(exevs)) # categs b.cuauhmaitl
Ftest.global <- matrix(NA, length(exevs), 5)# + chingon global
Ftest.partial <- as.list(1:length(exevs))# + chingon partial
# from plspm
Y.lvs <- pls$scores
IDM <- pls$model$IDM
blocks <- pls$model$blocks
modes <- pls$model$modes
scheme <- pls$model$scheme
scaled <- pls$model$scaled
tol <- pls$model$tol
iter <- pls$model$iter
endo <- rowSums(IDM)
endo[endo!=0] <- 1
n.endo <- sum(endo)
# size limite
if (size < 1) size.limit=ceiling(size*N) else size.limit=size
# for each exev
for (i in 1:length(exevs))
{
if (exevs[i] == 1) # ya se uso
{
Ftest.global[i,4] <- 1 # sta bien chafa
Ftest.partial[[i]] <- matrix(1,1,5)# sta bien chafa
list.cata[[i]] <- NA# NAs a.cuauhmaitl
list.catb[[i]] <- NA# NAs b.cuauhmaitl
next
} else # se usa
{
### moxexeloa
cat.exe <- unique(E.node[,i])
if (type.exev[i] == "ord") {
bin.split <- .ordinal.split(sort(cat.exe))
} else {# nom
if (exevs[i]==2) {
bin.split <- as.list(cat.exe)
} else {
bin.split <- .nominal.split(cat.exe)
}
}
### sanchez-aluja moxexeloa
Ftest.glo <- matrix(NA, length(bin.split[[1]]), 5) # global
Ftest.par <- as.list(1:length(bin.split[[1]])) # partial
for (aux in 1:length(bin.split[[1]])) # p/c moxexeloa
{
split.a <- which(E.node[,i] %in% bin.split[[1]][[aux]])# a.cuauhmaitl
split.b <- which(E.node[,i] %in% bin.split[[2]][[aux]]) # b.cuauhmaitl
n.a <- length(split.a)
n.b <- length(split.b)
if (n.a<=size.limit || n.b<=size.limit)# cuahuitl aturat
{
Ftest.glo[aux,4] <- 1 # sta bien chafa
Ftest.par[[aux]] <- matrix(1,1,5)# sta bien chafa
next
} else# cuahuitl segueix
{
idm <- .sanchez.aluja(Y.lvs, IDM, split.a, split.b)
Ftest.glo[aux,] <- idm[[1]] # global
Ftest.par[[aux]] <- idm[[2]] # partial
}
}
### el gallo dels exevs
if (mox=="pathmox")
min.p <- which(Ftest.glo[,4]==min(Ftest.glo[,4]))[1]
if (mox=="techmox")
{
pvals.ave <- unlist(lapply(Ftest.par, function(x) exp(mean(log(x[,4]))) ))
min.p <- which(pvals.ave==min(pvals.ave))[1]
}
Ftest.global[i,] <- Ftest.glo[min.p,]
Ftest.partial[[i]] <- Ftest.par[[min.p]]
list.cata[[i]] <- bin.split[[1]][[min.p]]
list.catb[[i]] <- bin.split[[2]][[min.p]]
}
}
### els chafillas
if (mox=="pathmox")
{
otras <- order(Ftest.global[,4])
F.otras <- Ftest.global[otras,]
F.otras[,4] <- round(F.otras[,4],6)
F.otras[,5] <- round(F.otras[,5],4)
colnames(F.otras) <- c("f.stat","df.num","df.den","p.val","sa/sb")
cata.otras <- NULL
catb.otras <- NULL
for (i in 1:length(otras))
{
cata.otras <- rbind(cata.otras, paste(as.vector(list.cata[[otras[i]]]),sep="",collapse="/"))
catb.otras <- rbind(catb.otras, paste(as.vector(list.catb[[otras[i]]]),sep="",collapse="/"))
}
test.otras <- data.frame(variable=as.vector(colnames(EXEV)[otras]), F.otras,
categ.a=cata.otras, categ.b=catb.otras)
test.otras <- test.otras[which(!is.na(F.otras[,1])),]
test.otras[,5] <- format(test.otras[,5], scientific=FALSE)
colnames(test.otras)[6] <- "sa/sb"
}
if (mox=="techmox")
{
pvals.ave <- unlist(lapply(Ftest.partial, function(x) exp(mean(log(x[,4]))) ))
otras <- order(pvals.ave)
cata.otras <- NULL
catb.otras <- NULL
for (i in 1:length(otras))
{
cata.otras <- rbind(cata.otras, paste(as.vector(list.cata[[otras[i]]]),sep="",collapse="/"))
catb.otras <- rbind(catb.otras, paste(as.vector(list.catb[[otras[i]]]),sep="",collapse="/"))
}
test.otras <- data.frame(variable=as.vector(colnames(EXEV)[otras]), pval.geomean=sort(pvals.ave),
categ.a=cata.otras, categ.b=catb.otras)
test.otras <- test.otras[which(!is.na(Ftest.global[,1])),]
test.otras[,2] <- format(test.otras[,2], scientific=FALSE)
}
### el chingon
pvals.ave <- unlist(lapply(Ftest.partial, function(x) exp(mean(log(x[,4]))) ))
if (mox=="pathmox")
optim <- which(Ftest.global[,4]==min(Ftest.global[,4]))[1]
if (mox=="techmox")
optim <- which(pvals.ave==min(pvals.ave))[1]
f.optim <- Ftest.global[optim,]
cats.optim <- list.cata[[optim]]
anticat.optim <- list.catb[[optim]]
aux.a <- which(E.node[,optim] %in% cats.optim)
aux.b <- which(E.node[,optim] %in% anticat.optim)
part.a <- indelem[aux.a,2]
part.b <- indelem[aux.b,2]
list.elems <- list(part.a, part.b)
### resul
test.global <- data.frame(exev=colnames(EXEV)[optim], f.stat=f.optim[1], df.num=f.optim[2],
df.den=f.optim[3], p.val=f.optim[4], s1_s2=round(f.optim[5],4))
colnames(test.global)[6] <- "sa/sb"
test.partial <- list(exev=colnames(EXEV)[optim], p.val=pvals.ave[optim],
F.test=Ftest.partial[[optim]])
categs <- list(as.vector(cats.optim), as.vector(anticat.optim))
resul <- list(inner.global=test.global, inner.partial=test.partial, categs=categs,
otras=test.otras, list.elems=list.elems)
return(resul)
}
.fix.xexeloa <-
function(pls, DT, EXEV, type.exev, elemnod, nv, size, mox)
{
# cuauhmaitl "nv" elems
N <- nrow(EXEV)
if (nv==0) elems<-1:length(elemnod) else elems<-which(elemnod==nv)
indelem <- cbind(1:length(elems), elems)# identificador
E.node <- EXEV[elems,]# elems en cuauhmaitl actual
exevs <- apply(E.node, 2, function(x) nlevels(as.factor(x)))# vector fo exevs
list.cata <- as.list(1:length(exevs)) # categs a.cuauhmaitl
list.catb <- as.list(1:length(exevs)) # categs b.cuauhmaitl
Ftest.global <- matrix(NA, length(exevs), 5)# + chingon global
Ftest.partial <- as.list(1:length(exevs))# + chingon partial
# parameters of plspm
Y.lvs <- pls$scores# matrix of latent variables
IDM <- pls$model$IDM# IDM matrix
blocks <- pls$model$blocks
modes <- pls$model$modes
scheme <- pls$model$scheme
scaled <- pls$model$scaled
tol <- pls$model$tol
iter <- pls$model$iter
endo <- rowSums(IDM)
endo[endo!=0] <- 1
n.endo <- sum(endo)
# size limite
if (size < 1) size.limit=ceiling(size*N) else size.limit=size
# for each exev
for (i in 1:length(exevs))
{
if (exevs[i] == 1) # ya se uso
{
Ftest.global[i,4] <- 1 # sta bien chafa
list.cata[[i]] <- NA# NAs a.cuauhmaitl
list.catb[[i]] <- NA# NAs b.cuauhmaitl
next # next exev
} else # se usa
{
### moxexeloa
cat.exe <- unique(E.node[,i])
if (type.exev[i] == "ord") {
bin.split <- .ordinal.split(sort(cat.exe))
} else { # nom
if (exevs[i]==2) {
bin.split <- as.list(cat.exe)
} else {
bin.split <- .nominal.split(cat.exe) # function of nominal splits
}
}
### sanchez-aluja moxexeloa
Ftest.glo <- matrix(NA, length(bin.split[[1]]), 5) # global
Ftest.par <- as.list(1:length(bin.split[[1]])) # partial
for (aux in 1:length(bin.split[[1]])) # p/c moxexeloa
{
split.a <- which(E.node[,i] %in% bin.split[[1]][[aux]])# a.cuauhmaitl
split.b <- which(E.node[,i] %in% bin.split[[2]][[aux]]) # b.cuauhmaitl
n.a <- length(split.a)
n.b <- length(split.b)
if (n.a<=size.limit || n.b<=size.limit)# cuahuitl aturat
{
Ftest.glo[aux,4] <- 1 # sta bien chafa
Ftest.par[[aux]] <- matrix(1,1,5)# sta bien chafa
next
} else# cuahuitl segueix
{
idm <- .sanchez.aluja(Y.lvs, IDM, split.a, split.b)
Ftest.glo[aux,] <- idm[[1]] # global
Ftest.par[[aux]] <- idm[[2]] # partial
}
}
### el gallo dels exev
if (mox=="pathmox")
min.p <- which(Ftest.glo[,4]==min(Ftest.glo[,4]))[1]
if (mox=="techmox")
{
pvals.ave <- unlist(lapply(Ftest.par, function(x) exp(mean(log(x[,4]))) ))
min.p <- which(pvals.ave==min(pvals.ave))[1]
}
Ftest.global[i,] <- Ftest.glo[min.p,]
Ftest.partial[[i]] <- Ftest.par[[min.p]]
list.cata[[i]] <- bin.split[[1]][[min.p]]
list.catb[[i]] <- bin.split[[2]][[min.p]]
}
}
### els chafillas
if (mox=="pathmox")
{
otras <- order(Ftest.global[,4])
F.otras <- Ftest.global[otras,]
F.otras[,4] <- round(F.otras[,4],6)
F.otras[,5] <- round(F.otras[,5],4)
colnames(F.otras) <- c("f.stat","df.num","df.den","p.val","sa/sb")
cata.otras <- NULL
catb.otras <- NULL
for (i in 1:length(otras))
{
cata.otras <- rbind(cata.otras, paste(as.vector(list.cata[[otras[i]]]),sep="",collapse="/"))
catb.otras <- rbind(catb.otras, paste(as.vector(list.catb[[otras[i]]]),sep="",collapse="/"))
}
test.otras <- data.frame(variable=as.vector(colnames(EXEV)[otras]), F.otras,
categ.a=cata.otras, categ.b=catb.otras)
test.otras <- test.otras[which(!is.na(F.otras[,1])),]
test.otras[,5] <- format(test.otras[,5], scientific=FALSE)
colnames(test.otras)[6] <- "sa/sb"
who.otras <- otras[which(!is.na(F.otras[,1]))]
}
if (mox=="techmox")
{
pvals.ave <- unlist(lapply(Ftest.partial, function(x) exp(mean(log(x[,4]))) ))
otras <- order(pvals.ave)
pvals.otras <- pvals.ave[otras]
cata.otras <- NULL
catb.otras <- NULL
for (i in 1:length(otras))
{
cata.otras <- rbind(cata.otras, paste(as.vector(list.cata[[otras[i]]]),sep="",collapse="/"))
catb.otras <- rbind(catb.otras, paste(as.vector(list.catb[[otras[i]]]),sep="",collapse="/"))
}
test.otras <- data.frame(variable=as.vector(colnames(EXEV)[otras]), pval.geomean=sort(pvals.ave),
categ.a=cata.otras, categ.b=catb.otras)
test.otras <- test.otras[which(!is.na(Ftest.global[,1])),]
test.otras[,2] <- format(test.otras[,2], scientific=FALSE)
who.otras <- otras[which(pvals.otras!=1)]
}
### el fijado
if (length(who.otras)>0)
{
candi.table <- cbind(Num=who.otras, test.otras)
if (nv==0)
cat("CANDIDATE SPLITS: ROOT NODE", "\n\n")
if (nv>0)
cat("CANDIDATE SPLITS: NODE", nv, "\n\n")
print(candi.table)
cat("\n")
cat("Select a candidate split from column 'Num' and then press Enter:", "\n")
fijo <- scan(what=integer(1), n=1)
if (mode(fijo)!="numeric" || length(fijo)!=1 || fijo<1 || (fijo%%1)!=0) {
cat("Invalid number of candidate split", "\n")
stop()
}
f.fijo <- Ftest.global[fijo,]
} else {
fijo <- which(Ftest.global[,4]==min(Ftest.global[,4]))[1]
f.fijo <- Ftest.global[fijo,]
}
cats.fijo <- list.cata[[fijo]]
anticat.fijo <- list.catb[[fijo]]
aux.a <- which(E.node[,fijo] %in% cats.fijo)
aux.b <- which(E.node[,fijo] %in% anticat.fijo)
part.a <- indelem[aux.a,2]
part.b <- indelem[aux.b,2]
list.elems <- list(part.a, part.b)
### resul
if (mox=="pathmox")
{
test.fijo <- data.frame(exev=colnames(EXEV)[fijo], f.stat=f.fijo[1], df.num=f.fijo[2],
df.den=f.fijo[3], p.val=f.fijo[4], sa_sb=round(f.fijo[5],4))
colnames(test.fijo)[6] <- "sa/sb"
}
if (mox=="techmox")
{
test.fijo <- list(exev=colnames(EXEV)[fijo], p.val=pvals.ave[fijo],
F.test=Ftest.partial[[fijo]])
}
categs <- list(as.vector(cats.fijo), as.vector(anticat.fijo))
resul <- list(inner.test=test.fijo, categs=categs, list.elems=list.elems, otras=test.otras)
return(resul)
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.