我正在处理分段包,在函数中调用davies.test()时遇到了问题。
考虑以下情况:
library(segmented)
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data = data)
fit.seg = segmented(fit, seg.Z = ~ x)
davies.test(fit.seg, seg.Z = ~ x, alternative = "greater")这是完美的工作,并表明分段回归有两个统计上不同的斜率。
现在,如果我将所有这些都打包成这样的函数:
testit <- function() {
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data)
fit.seg = segmented(fit, seg.Z = ~ x)
davies.test(fit.seg, seg.Z = ~ x, alternative = "greater")$p.value
}
testit()那就行了..。
但是如果我从全局环境中删除fit,那么它就失败了。
> rm(fit)
> testit()
Error in eval(expr, envir, enclos) : object 'fit' not found问题似乎在于davies.test试图访问封装在fit中的数据的方式:它似乎没有在封闭范围(在本例中是testit函数)中寻找fit,而是直接跳到全局范围。
我确信这个问题与R的范围规则的一些微妙之处有关。如果我能找到一个可以防止我用这种边缘情况困扰包作者的快速修复方法,那就太好了。
谢谢安德鲁。
发布于 2016-01-04 13:46:50
尝试插入下面标有##的行。仍然存在这样的差异,正如运行修改后的testit时出现的警告所示,但是输出pvalue是相同的,因此它可能足以满足您的需要。当然,这是包中的一个bug,最好是询问包的维护人员是否会修复它。
library(segmented)
testit <- function() {
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data)
fit.seg = segmented(fit, seg.Z = ~ x)
environment(davies.test) <- environment() ##
davies.test(fit.seg, seg.Z = ~ x, alternative = "greater")$p.value
}
testit()给予:
[1] 0.01858149
Warning message:
In summary.lm(object) : essentially perfect fit: summary may be unreliable发布于 2016-01-04 13:40:19
不需要让它成为一个全局变量。问题实际上是在segmented,而不是davies.test。这不是找到fit。
可以使用dynGet在任何环境(包括调用函数的环境)中定位fit:
testit <- function() {
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data)
fit.seg = segmented(dynGet("fit"), seg.Z = ~ x)
davies.test(fit.seg, seg.Z = ~ x, alternative = "greater")$p.value
}
testit()这应该能按你的意愿运作。
如果在不同的环境中有多个名为fit的变量,那么使用get (参见?get)来指定要从哪个环境获取它。dynGet是“随处可见,先返回”的懒散版本。
发布于 2016-01-07 05:57:46
我联系了segmented的作者,他立即回复了。他对原问题提出的另一个解决办法是
testit <- function() {
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data)
fit.seg = segmented(fit, seg.Z = ~ x)
fit.seg$call$obj<-fit
davies.test(fit.seg, seg.Z = ~ x, alternative = "greater")$p.value
}但是,他还指出,lm对象实际上应该直接传递给davies.test(),如下所示:
testit <- function() {
data = data.frame(x = 1:21, y = c(10:1, 0:10))
fit = lm(y ~ x, data)
davies.test(fit, seg.Z = ~ x, alternative = "greater")$p.value
}但是,为了澄清这一点,应该注意到这两段代码做了不同的事情:第二个片段实际上实现了我的最初目的(检查fit中是否有统计意义的中断),而第一个片段检查是否有第二个中断。
https://stackoverflow.com/questions/34591647
复制相似问题