大家好,在S中偏相关系数用什么函数来算啊?

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

大家好,有人知道在SPLUS中偏相关系数怎么算啊?谢谢啦

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏0 回帖

关键词：偏相关系数 相关系数 偏相关 大家好 PLUS 函数 系数 偏相关 皮尔逊相关系数 偏相关系数 pearson相关系数 相关系数检验 spearman相关系数 复相关系数 相关系数矩阵

> x=rnorm(100)
> x.pacf=pacf(x)
> x.pacf

Partial autocorrelations of series 'x', by lag

1 2 3 4 5 6 7 8 9 10 11
-0.062 0.009 0.042 -0.290 -0.155 0.017 -0.010 -0.180 0.001 -0.010 -0.044
12 13 14 15 16 17 18 19 20
-0.089 -0.026 0.070 0.066 -0.049 -0.096 0.103 -0.085 0.070

pacf()函数在R中才有的，在S中没有这个函数啊

> pacf(x)
Problem: Couldn't find a function definition for "pacf"
Use traceback() to see the call stack
我的意思是控制了X3的情况下计算X1和X2的相关系数啊。

该怎么写呢？谢谢啦

不好意思,是我的理解错误,将你的偏相关系数理解成了偏自相关系数.

这是两个不同的概念.

发现这个问题之后, 我也在想如何计算偏相关系数, 找了好的资料,呵呵都没有讲这个问题.

发现SPLUS 或R 中没有现成的函数可以计算它们.

但是还是有一些人自己写code来计算. 在R中

John Fox 的代码:

Suppose that R is the correlation matrix among the variables. Then

Rinv <- solve(R)
D <- diag(1 / sqrt(diag(Rinv)))
P <- -D %*% Rinv %*% D

Martyn Plummer:

In general you invert the variance-covariance matrix and then rescale it so the diagonal is one. The off-diagonal elements are the negative partial correlation coefficients given all other variables.

pcor2 <- function(x){
conc <- solve(var(x))
resid.sd <- 1 / sqrt(diag(conc))
pcc <- -sweep(sweep(conc, 1, resid.sd, "*"), 2, resid.sd, "*")
return(pcc)
}
pcor2(cbind(x1, x2, x3))

KP responded: This is the version I’m using now, together with a test for significance of each coefficient (H0: coeff = 0):

f.parcor <- function (x, test = FALSE, p = 0.05) {
nvar <- ncol(x)
ndata <- nrow(x)
conc <- solve(cor(x))
resid.sd <- 1 / sqrt(diag(conc))
pcc <- -sweep(sweep(conc, 1, resid.sd, "*"), 2, resid.sd, "*")
colnames(pcc) <- rownames(pcc) <- colnames(x)
if (test) {

t.df <- ndata - nvar
t <- pcc / sqrt((1 - pcc^2) / t.df)
pcc <- list(coefs = pcc, significant = t > qt(1 - (p/2),
df = t.df))
}
return(pcc)
}

[此贴子已经被作者于2007-4-14 9:36:10编辑过]

好的,谢谢,让我好好看一下代码先