3
$\begingroup$

I am looking for a R command to test the difference of two linear regressoon betas. Lets say I have data $x_1, x_2...x_{n+1}$.

$\beta_1$ is obtained from regressing $x_1$ to $x_n$ onto $1$ to $n$.

$\beta_2$ is obtained from regressing $x_2$ to $x_{n+1}$ onto $1$ to $n$.

Is there a way in R to test whether $\beta_1$ and $\beta_2$ are statistically different?

Edit:

Here is a clearer problem statement:

I changed the notation for data from $x$ to $z$...

That's it. Should be very clear now... Thanks!

Data observations: $z_1, z_2, ..., z_{n+1}$

y1 = z_1,z_2,.........  z_n 
y2 = z_2, z_3,......... z_{n+1}

x1 = 1, ..., n
x2 = 1, ..., n

y1 = A1+ x1 * B1 + epsilon_1
y2 = A2 + x2 * B2 + epsilon_2

H0: B1 and B2 are statistically significally different...

$\endgroup$

closed as off-topic by Nick Cox, gung, Peter Flom Sep 20 '16 at 21:26

This question appears to be off-topic. The users who voted to close gave this specific reason:

  • "This question appears to be off-topic because EITHER it is not about statistics, machine learning, data analysis, data mining, or data visualization, OR it focuses on programming, debugging, or performing routine operations within a statistical computing platform. If the latter, you could try the support links we maintain." – Nick Cox, gung, Peter Flom
If this question can be reworded to fit the rules in the help center, please edit the question.

  • $\begingroup$ I'm a bit confused by your notation. Is "x1, x2...x(n+1)" referring to random variables or observations? And what is the significance of "(n+1)"? $\endgroup$ – Andy W Jan 27 '12 at 16:02
  • $\begingroup$ Presumably you mean "H0: B1 and B2 are the same" and you are looking for statistically significant evidence to reject H0. $\endgroup$ – Peter Ellis Jan 30 '12 at 23:58
  • $\begingroup$ Are you trying to model a structural change in the data? If so, then the package strucchange might be appropriate. $\endgroup$ – Sean Jan 31 '12 at 9:04
4
$\begingroup$

You can stack the data so that you have $x_1 ... x_n$ followed by $x_2 ... x_{n+1}$ in one column, then the next column would have $1 ... n$ repeated twice, then a third column would have $0$'s coresponding to the 1st group and $1$'s for the second group. Now do the regression including an interaction between the $1 ... n$ and the $0/1$ group variable. The slope for the interaction represents the difference between the 2 slopes you are interested in, testing the interaction will test if the slopes are different. However, the standard regression assumptions may not hold here (your data certainly is not independent), so you should do some type of simulation under the null hypothesis to determine the critical region to use to determine significance.

Here is some example R code:

y1 <- anscombe$y1[order(anscombe$x1)]
y2 <- anscombe$y2[order(anscombe$x2)]


df1 <- data.frame( y=c(y1[-11], y1[-1]), x=rep(1:10, 2), g=rep(0:1, each=10))
df2 <- data.frame( y=c(y2[-11], y2[-1]), x=rep(1:10, 2), g=rep(0:1, each=10))

fit1 <- lm( y ~ x*g, data=df1 )
fit2 <- lm( y ~ x*g, data=df2 )

tmpfun <- function(n, beta, sigma) {
    x <- 1:n
    y <- beta*x + rnorm(n,0,sigma)
    df <- data.frame( y=c(y[-n],y[-1]), x=rep(seq(to=n-1), 2), 
        g=rep(0:1, each=n-1) )
    fit <- lm( y~x*g, data=df )
    coef(fit)[4]
}

tmpfit <- lm(y1 ~ seq_along(y1))
out1 <- replicate(1000, tmpfun(11, coef(tmpfit)[2], summary(tmpfit)$sigma) )
hist(out1)
abline( v=quantile(out1, c(0.025, 0.975)), col='red')
abline( v=coef(fit1)[4], col='blue' )

tmpfit <- lm(y2 ~ seq_along(y1))
out2 <- replicate(1000, tmpfun(11, coef(tmpfit)[2], summary(tmpfit)$sigma) )
hist(out2)
abline( v=quantile(out2, c(0.025, 0.975)), col='red')
abline( v=coef(fit2)[4], col='blue' )

If you tell us more about what you are trying to accomplish we may be able to suggest better methods.

$\endgroup$
  • $\begingroup$ Sorry for not being clear in my previous problem statement: I changed the notation for data from x to z... That's it. Should be very clear now... Thanks! Data: z1, z2, ..., z_{n+1} y1 = z_1,z_2,......... z_n y2 = z_2, z_3,......... z_{n+1} x1 = 1, ..., n x2 = 1, ..., n y1 = A1+ x1 * B1 + epsilon_1 y2 = A2 + x2 * B2 + epsilon_2 H0: B1 and B2 are statistically significally different... $\endgroup$ – Luna Jan 27 '12 at 21:11
  • 1
    $\begingroup$ "Stacking" the data in this manner is clever. However, that induces perfect correlations among many of the responses. Because your model does not appear to account for that, its estimates look problematic--and even if they are correct, the covariance matrix, p-values, intervals, and associated hypothesis tests will certainly be wrong. This problem is likely so severe that "your data certainly is not independent" doesn't look like a sufficiently strong caution. $\endgroup$ – whuber Sep 20 '16 at 18:53
  • $\begingroup$ @Greg Snow: Can you explain what out1 & out2 in your code are? $\endgroup$ – user603 Oct 11 '17 at 22:00
  • $\begingroup$ @user603, The tmpfun function returns the coefficient for the interaction in the analysis of the simulated data. out1 and out2 are vectors with the returns from tmpfun from 1,000 runs of the simulation (so vectors of interaction coefficients from 1,000 simulation runs). $\endgroup$ – Greg Snow Oct 12 '17 at 13:41
  • $\begingroup$ @GregSnow; thanks. What is the difference between them? I mean why do you have two? $\endgroup$ – user603 Oct 12 '17 at 15:13
1
$\begingroup$

Are you doing something like:

set.seed(1)

x1<-rnorm(100,1) 
x2<-rnorm(100,mean(x1))
y<-rnorm(100, mean(x1+x2)+rnorm(100))

fit1<-lm(y~x1, data)
fit2<-lm(x1~x2, data)

If so, you might be able to get away with:

d<-coef(fit1)[2] - coef(fit2)[2] 

coef(fit)[2] is the second coefficient in the model, i.e. the beta

var1<-summary(fit1)$coef[4]^2 

This is the fourth coefficient in the summary, or the se of beta1

var2<-summary(fit2)$coef[4]^2


ztest<- d / sqrt(var1+var2)

1-pnorm(ztest)

       x1 
0.5918211 

this was my answer, so I would say there is no difference, but beta1 was -.11 and beta 2 was .-07, so they're pretty similar.

$\endgroup$
  • 1
    $\begingroup$ Welcome to the site Corey. Please use the markdown formatting available (see the buttons at the top of editor when inserting your answer), or take a look at the advanced markdown editing help, stats.stackexchange.com/editing-help. In particular, it is much easier to read code when it is formatted appropriately. $\endgroup$ – Andy W Jan 27 '12 at 16:37

Not the answer you're looking for? Browse other questions tagged or ask your own question.