On 2011-03-11 17:10, jonbfish wrote:
Thanks again. I guess I have been fitting in splus gui and in R Commander, not
really calling it myself. I was thinking maybe there was a way to define it as
a class or something.
I will have to look at segments. Any thoughts if this might be easy to use with
years on an axis?
Alright, I see what you're after. (It would have helped
if you had mentioned R Commander sooner.)
All you need is to use lines() with appropriate endpoints.
Here's how:
x <- 1999:2010
y <- c(17.942,3.43,14.062,14.814,13.778,13.706,
9.748,13.088,12.1309728,9.644646671,9.134,8.84)
mod <- lm(y ~ x)
kt <- function(x, y){
yy <- outer(y, y, "-")
xx <- outer(x, x, "-")
z <- yy / xx
s <- z[lower.tri(z)]
slope <- median(s)
intercept <- median(y) - slope * median(x)
list(intercept = intercept, slope = slope)
}
ind1 <- which.min(x)
ind2 <- which.max(x)
w <- kt(x, y)
a <- w[["intercept"]]
b <- w[["slope"]]
x.ends <- c(x[ind1], x[ind2])
y.ends.kt <- a + b * x.ends
a <- coef(mod)[1]
b <- coef(mod)[2]
y.ends.lm <- a + b * x.ends
plot(y ~ x)
lines(x.ends, y.ends.kt, col = "blue")
lines(x.ends, y.ends.lm, col = "red")
Peter Ehlers
From: Peter Ehlers [via R]
[mailto:ml-node+3349445-1733968405-216...@n4.nabble.com]
Sent: Friday, March 11, 2011 06:38 PM
To: Anthony Seeman
Subject: Re: Kendall Theil line as fit?
On 2011-03-11 14:43, jonbfish wrote:
Thanks for the response, sorry I didn't post it initially.
kt.mat<-
function(x,y,z){
for(i in 1:length(x)){for(j in
1:length(y)){z[i,j]<-(y[j]-y[i])/(x[j]-x[i])}}
return(z)}
kt.slope<-
function(x,y,z,s){
count<-0
for(i in 1:length(x)){for(j in 1:length(y)){
if(j>= i+1) {
count<-count+1
s[count]<-z[i,j]}
}}
print(count)
return(s)}
#Site23
x<- c(1999,2000,2001,2002,2003,2004,2005,2006,2007,2008,2009,2010)
y<-
c(17.942,3.43,14.062,14.814,13.778,13.706,9.748,13.088,12.1309728,9.644646671,9.134,8.84)
z<-matrix(0:0,length(x),length(y))
z<-kt.mat(x,y,z)
z
s<-c(1:(length(x)*(length(x)-1)/2))
s<-kt.slope(x,y,z,s)
s
slope=median(s)
intercept=median(y)-slope*median(x)
cbind(slope,intercept)
plot(x,y)
abline(intercept,slope)
Okay, you're using abline() for the KT line.
But I still don't know what you're after.
From your original post:
Is there a way to make it appear like a regression fit instead
of a line that extends from the edges of the plot? I would like
to have the OLS appear as a dotted line and the KT a solid line
but as it is the KT line is longer.
So how are plotting your 'regression fit'?
abline( lm( y ~ x ) ) would also extend across the plot.
I suppose that you could use segments() with the
range of x-values.
BTW, here's a shorter version of your code:
yy<- outer(y, y, "-")
xx<- outer(x, x, "-")
z<- yy / xx
s<- z[lower.tri(z)]
slope<- median(s)
intercept<- median(y) - slope * median(x)
Peter Ehlers
______________________________________________
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.
