我想模仿以下的一段R代码统计我的数据(名为“mydata”)
在输入语言后,系统提示错误(红色字样)。由于是初学,不知问题出在哪里。本来想高分悬赏,可是我连一个论坛币也没有。因毕业在即,十分焦虑。请哪位好心人指点一下。万分感激!
library(pls)library(GeneNet)library(graph)library(Rgraphviz) example = function() {
data(ecoli)
# Apply PLS-PC network using 3 components on ecoli data
# Plot the resulting graph
result = PLSPC.network(ecoli)
} PLSPC.network = function(X, nbre = 3, method = "kernelpls", viz = 1){# X : Microarray experiments (n * p)#
n = number of experiments#
p = number of genes# nbre = Number of PLS components to include in the PLSPC network (default = 3 components)# method = implementation of the PLS algorithm (default = Kernel PLS)# viz = visualization (default = 1 -> visualization desired)
# Output : list of variables #
pcor
: PLSPC matrix #
graph : data.frame with the significant edges#
gr
: a graph object#
coef
: PLS matrix coefficients if (viz ==1){X = X[ , as.logical(1-duplicated(colnames(X)))]X = X[ , which(colnames(X) != "NA")]} if ((nbre > qr(scale(X, center = TRUE, scale = TRUE))$rank) && (is.null(nbre) == FALSE)){ print(paste("Number of PLS components must be < ", qr(scale(X, center = TRUE, scale = TRUE))$rank))} else{# Construction of the Partial Correlation matrix
Package PLS required
# pcor = coef = diag(ncol(X))for (i in 1:ncol(X)){
temp = mvr(X[, i]~X[, -i], ncomp = nbre, method = method, scale = TRUE, stripped = TRUE)
coef[i, -i] = as.data.frame(temp$coefficients)[ , nbre]} for (i in 1:ncol(X)){
for (j in i:ncol(X)){
if (i==j) pcor[i, j] = 1
else {
if(sign(coef[j, i]) != sign(coef[i, j])) pcor[i, j] = 0
else{
pcor[i, j] = sign(coef[j, i])*sqrt(coef[i, j]*coef[j, i])
}
}
pcor[j, i] = pcor[i, j]
}
}}if (viz == 1){# Local fdr multiple testing on the PLS
test.results = ggm.test.edges(pcor, plot = FALSE)signif = test.results$prob > 0.8
if (sum(signif)>1){# Visualization of the PLSPC network node.labels = colnames(X)gr = ggm.make.graph(test.results[signif, ], node.labels, drop.singles = FALSE)x11(); plot(gr,"twopi")return(list(pcor = pcor, graph = test.results[signif, ], gr = gr, coef = coef))}else{print("No significant interaction !");return(list(pcor = pcor, coef = coef))}}else{return(list(pcor = pcor, coef = coef))}}我载入程序包并按照代码输入语言后的界面显示如下:
> example = function()
+ {
data(ecoli)
+
result = PLSPC.network(ecoli)
> PLSPC.network = function(X, nbre = 3, method = "kernelpls", viz = 1){
+ if (viz ==1){
+ X = X[ , as.logical(1-duplicated(colnames(X)))]
+ X = X[ , which(colnames(X) != "NA")]
+ }
+
+ if ((nbre > qr(scale(X, center = TRUE, scale = TRUE))$rank) && (is.null(nbre) == FALSE)){
+ print(paste("Number of PLS components must be < ", qr(scale(X, center = TRUE, scale = TRUE))$rank))
+ }
+ else{
+ pcor = coef = diag(ncol(X))
+ for (i in 1:ncol(X)){
+ temp = mvr(X[, i]~X[, -i], ncomp = nbre, method = method, scale = TRUE, stripped = TRUE)
+
coef[i, -i] = as.data.frame(temp$coefficients)[ , nbre]
+ }
+ for (i in 1:ncol(X)){
+
for (j in i:ncol(X)){
+
if (i==j) pcor[i, j] = 1
+
else {
+
if(sign(coef[j, i]) != sign(coef[i, j])) pcor[i, j] = 0
+
else{
+
pcor[i, j] = sign(coef[j, i])*sqrt(coef[i, j]*coef[j, i])
+
}
+
}
+
pcor[j, i] = pcor[i, j]
+
}
+
}
+ }
+ if (viz == 1){
+ test.results = ggm.test.edges(pcor, plot = FALSE)
+ signif = test.results$prob > 0.8
+ if (sum(signif)>1){
+ node.labels = colnames(X)
+ gr = ggm.make.graph(test.results[signif, ], node.labels, drop.singles = FALSE)
+ x11(); plot(gr,"twopi")
+ return(list(pcor = pcor, graph = test.results[signif, ], gr = gr, coef = coef))}
+ else{print("No significant interaction !");return(list(pcor = pcor, coef = coef))}
+ }
+ else{return(list(pcor = pcor, coef = coef))}
+ }
> d=read.table("mydata.txt")
> result = PLSPC.network(d) (此两行为我与代码不同之处)
错误于model.frame.default(formula = X[, i] ~ X[, -i]) :
变元'X[, -i]'的种类(list)不对