分段程序包中的错误:断点混乱 [英] Errors in segmented package: breakpoints confusion

查看:65
本文介绍了分段程序包中的错误:断点混乱的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

使用分段包创建分段线性回归时,尝试设置自己的断点时出现错误;似乎只有当我尝试设置两个以上时.

Using the segmented package to create a piecewise linear regression I am seeing an error when I try to set my own breakpoints; it seems only when I try to set more than two.

(编辑)这是我正在使用的代码:

(EDIT) Here is the code I am using:

# data
bullard <- structure(list(Rt = c(0, 4.0054, 25.1858, 27.9998, 35.7259, 39.0769, 
45.1805, 45.6717, 48.3419, 51.5661, 64.1578, 66.828, 111.1613, 
114.2518, 121.8681, 146.0591, 148.8134, 164.6219, 176.522, 177.9578, 
180.8773, 187.1846, 210.5131, 211.483, 230.2598, 262.3549, 266.2318, 
303.3181, 329.4067, 335.0262, 337.8323, 343.1142, 352.2322, 367.8386, 
380.09, 388.5412, 390.4162, 395.6409), Tem = c(15.248, 15.4523, 
16.0761, 16.2013, 16.5914, 16.8777, 17.3545, 17.3877, 17.5307, 
17.7079, 18.4177, 18.575, 19.8261, 19.9731, 20.4074, 21.2622, 
21.4117, 22.1776, 23.4835, 23.6738, 23.9973, 24.4976, 25.7585, 
26.0231, 28.5495, 30.8602, 31.3067, 37.3183, 39.2858, 39.4731, 
39.6756, 39.9271, 40.6634, 42.3641, 43.9158, 44.1891, 44.3563, 
44.5837)), .Names = c("Rt", "Tem"), class = "data.frame", row.names = c(NA, 
-38L))

library(segmented)

# create a linear model
out.lm <- lm(Tem ~ Rt, data=bullard)

o<-segmented(out.lm, seg.Z=~Rt, psi=list(Rt=c(200,300)), control=seg.control(display=FALSE))

使用psi选项,我尝试了以下操作:

Using the psi option, I have tried the following:

psi = list(x = c(150, 300)) -- OK
psi = list(x = c(100, 200)) -- OK
psi = list(x = c(200, 300)) -- OK
psi = list(x = c(100, 300)) -- OK
psi = list(x = c(120, 150, 300)) -- error 1 below
psi = list(x = c(120, 300)) -- OK
psi = list(x = c(120, 150)) -- OK
psi = list(x = c(150, 300)) -- OK
psi = list(x = c(100, 200, 300)) -- error 2 below

(1)Error in segmented.lm(out.lm, seg.Z = ~Rt, psi = list(Rt = c(120, 150, : only 1 datum in an interval: breakpoint(s) at the boundary or too close

(1) Error in segmented.lm(out.lm, seg.Z = ~Rt, psi = list(Rt = c(120, 150, : only 1 datum in an interval: breakpoint(s) at the boundary or too close

(2)Error in diag(Cov[id, id]) : subscript out of bounds

我已经在此问题上列出了我的数据,但作为指导,x数据的限制大约为0 --400.

I have already listed my data at this question, but as a guide the limits on the x data are about 0--400.

与此相关的第二个问题是:我实际上如何使用此分段程序包修复断点?

A second question that pertains to this one is: how do I actually fix the breakpoints using this segmented package?

推荐答案

此处的问题似乎是在segmented程序包中捕获错误的能力很差.查看segmented.lm的代码可以进行一些调试.例如,在psi = list(x = c(100, 200, 300))的情况下,拟合增强线性模型,如下所示:

The issue here seems to be poor error trapping in the segmented package. Having a look at the code for segmented.lm allows a bit of debugging. For example, in the case of psi = list(x = c(100, 200, 300)), an augmented linear model is fitted as shown below:

lm(formula = Tem ~ Rt + U1.Rt + U2.Rt + U3.Rt + psi1.Rt + psi2.Rt + 
    psi3.Rt, data = mf)

Call:
lm(formula = Tem ~ Rt + U1.Rt + U2.Rt + U3.Rt + psi1.Rt + psi2.Rt + 
    psi3.Rt, data = mf)

Coefficients:
(Intercept)           Rt        U1.Rt        U2.Rt        U3.Rt      psi1.Rt        
   15.34303      0.04149      0.04591    742.74186   -742.74499      1.02252       
   psi2.Rt      psi3.Rt  
        NA           NA  

如您所见,拟合具有NA值,然后这些值会生成退化的方差-协方差矩阵(在代码中称为Cov).该功能不会对此进行检查,而是尝试从Cov中拉出对角线条目并失败,并显示错误消息.至少第一个错误虽然可能不太有用,但还是被函数本身捕获了,表明断点太近了.

As you can see, the fit has NA values which then result in a degenerate variance-covariance matrix (called Cov in the code). The function doesn't check for this and tries to pull out diagonal entries from Cov and fails with the error message shown. At least the first error, although perhaps not overly helpful, is caught by the function itself and suggests that the break-points are too close.

在函数中没有更好的错误捕获的情况下,我认为您可以做的就是采用反复试验的方法(并避免断点过于靠近).例如,psi = list(x = c(50, 200, 300))似乎可以正常工作.

In the absence of better error trapping in the function, I think that all you can do is adopt a trial and error approach (and avoid break points which are too close). For example, psi = list(x = c(50, 200, 300)) seems to work ok.

这篇关于分段程序包中的错误:断点混乱的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆