My data:
I have raw data points that form a logit style curve as if they were a time
series. Which is to say they form 3 distinct lines with 3 distinct slopes
in backwards z pattern. A certain class of my data looks essentially flat
to the eye with marginal oscillation. What is important to me is the x
value at which the state change is occurring, in other words, the break
point
Use of segmented():
Segmented does a very good job of capturing the breakpoints and fitting
three distinct slopes, i.e. linear models. It does not, however give me
Pr(>|t|) for the break point coefficients. I need to answer the question
H:0 Beta0=Beta with a certainty metric, i.e. probability statistic. This is
especially important for my, flat looking data class.
davies.test() question:
davies.test() only excepts lm() or glm() objects as input. If I run
segmented to find 1 breakpoint instead of 2, I get a totally bogus answer.
Without knowing the breakpoints, how is this test able to assess the proper
breakpiont? It appears to only give 1 best breakpoint, which is not
consistent with the breakpoints found by segmented(). Also, is K the number
of breakpoints or the number of iterations that it evaluates the breakpoint?
Thanks in advance.
lmfit<-glm(TotRad_KW~HRRPUA_kWm2,data=d1)
davies.test(lmfit,seg.Z=~HRRPUA_kWm2,k=1000,alternative="less",
beta0=0,dispersion=NULL)
Davies' test for a change in the slope
data: Model = gaussian , link = identity
formula = TotRad_KW ~ HRRPUA_kWm2
segmented variable = HRRPUA_kWm2
`Best' at = 561.205, n.points = 1000, p-value < 2.2e-16
alternative hypothesis: less
segments <- segmented(lmfit, seg.Z=~HRRPUA_kWm2,psi=c(475,550))
summary(segments)
***Regression Model with Segmented Relationship(s)***
Call:
segmented.glm(obj = lmfit, seg.Z = ~HRRPUA_kWm2, psi = c(475,
550))
Estimated Break-Point(s):
Est. St.Err
psi1.HRRP 430.2 4.087
psi2.HRRP 484.6 3.077
t value for the gap-variable(s) V: 0 0
Meaningful coefficients of the linear terms:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -38.6993 274.7666 -0.141 0.8891
HRRPUA_kWm2 1.4297 0.7472 1.914 0.0668 .
U1.HRRP 42.2884 4.7696 8.866 NA
U2.HRRP -40.5897 4.7123 -8.614 NA
---
Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
(Dispersion parameter for gaussian family taken to be 6934.706)
Null deviance: 70776718 on 31 degrees of freedom
Residual deviance: 180302 on 26 degrees of freedom
AIC: 377.19
Convergence attained in 2 iterations with relative change -1.662839e-14
--
Greg Cohn
Forestry Technician
USDA Forest Service Fire Lab
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
