Вывод о коэффициенте наклона в R

По умолчанию lm суммарный коэффициент наклона теста равен нулю. Мой вопрос очень простой. Я хочу знать, как проверить коэффициент наклона, равный ненулевому значению. Одним из подходов может быть использование confint, но это не дает p-значения. Мне тоже интересно, как сделать односторонний тест с lm.

ctl <- c(4.17,5.58,5.18,6.11,4.50,4.61,5.17,4.53,5.33,5.14)
trt <- c(4.81,4.17,4.41,3.59,5.87,3.83,6.03,4.89,4.32,4.69)
group <- gl(2,10,20, labels=c("Ctl","Trt"))
weight <- c(ctl, trt)
lm.D9 <- lm(weight ~ group)
summary(lm.D9)

Call:
lm(formula = weight ~ group)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.0710 -0.4938  0.0685  0.2462  1.3690 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   5.0320     0.2202  22.850 9.55e-15 ***
groupTrt     -0.3710     0.3114  -1.191    0.249    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.6964 on 18 degrees of freedom
Multiple R-squared: 0.07308,    Adjusted R-squared: 0.02158 
F-statistic: 1.419 on 1 and 18 DF,  p-value: 0.249 


confint(lm.D9)
              2.5 %    97.5 %
(Intercept)  4.56934 5.4946602
groupTrt    -1.02530 0.2833003

Спасибо за ваше время и усилия.


r stat lm
person MYaseen208    schedule 11.11.2011    source источник
comment
вы можете увидеть в выводе summary(lm.D9).   -  person kohske    schedule 11.11.2011
comment
@kohske: Спасибо за ваш комментарий. summary обеспечивает проверку коэффициентов, равных нулю.   -  person MYaseen208    schedule 11.11.2011

Ответы (6)


arrow_upward
2
arrow_downward

Используйте функцию linearHypothesis из пакета car. Например, вы можете проверить, равен ли коэффициент groupTrt -1, используя.

linearHypothesis(lm.D9, "groupTrt = -1")

Linear hypothesis test

Hypothesis:
groupTrt = - 1

Model 1: restricted model
Model 2: weight ~ group

  Res.Df     RSS Df Sum of Sq      F  Pr(>F)  
1     19 10.7075                              
2     18  8.7292  1    1.9782 4.0791 0.05856 .
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
person Ramnath    schedule 11.11.2011
comment
Большое спасибо за ваш хороший ответ. Интересно, как сделать односторонний тест. Я пробовал это с groupTrt >= -1, но это не сработало. - person MYaseen208; 11.11.2011
comment
Как насчет уменьшения вдвое значения p... ?? (Это был бы стандартный ответ на вопрос, как получить одностороннее значение p, хотя вы должны проверить комментарии ниже ответа @kohske) - person Ben Bolker; 12.11.2011

arrow_upward
7
arrow_downward

как говорит @power, вы можете сделать это своими руками. вот пример:

> est <- summary.lm(lm.D9)$coef[2, 1]
> se <- summary.lm(lm.D9)$coef[2, 2]
> df <- summary.lm(lm.D9)$df[2]
> 
> m <- 0
> 2 * abs(pt((est-m)/se, df))
[1] 0.2490232
> 
> m <- 0.2
> 2 * abs(pt((est-m)/se, df))
[1] 0.08332659

и вы можете выполнить односторонний тест, опуская 2*.

ОБНОВЛЕНИЯ

вот пример двусторонней и односторонней вероятности:

> m <- 0.2
> 
> # two-side probability
> 2 * abs(pt((est-m)/se, df))
[1] 0.08332659
> 
> # one-side, upper (i.e., greater than 0.2)
> pt((est-m)/se, df, lower.tail = FALSE)
[1] 0.9583367
> 
> # one-side, lower (i.e., less than 0.2)
> pt((est-m)/se, df, lower.tail = TRUE)
[1] 0.0416633

обратите внимание, что сумма верхней и нижней вероятностей равна ровно 1.

person kohske    schedule 11.11.2011
comment
Я не думаю, что опускание 2 является статистически правильным. Я думаю, вам нужно подставить 1,644854 = qnorm (.95), а затем смотреть только в направлении, указанном последней еще не сформулированной гипотезой. - person IRTFM; 11.11.2011
comment
Это выглядит лучше. Может, я раньше не понимал, что ты делаешь. Мне удобнее смотреть на среднее + 1,64 * se (или среднее - 1,64 * se в зависимости от конкретной гипотезы) по сравнению с 0. Мне кажется, что большинство людей, задающих этот вопрос, имеют предельный результат, который они просто пытаются протолкнуть. искусственный финиш 0,05. - person IRTFM; 11.11.2011
comment
Или используйте t.test для данных напрямую; см. мой ответ ниже. - person James; 11.11.2011

arrow_upward
1
arrow_downward

В пакете smatr есть функция slope.test(), с помощью которой можно использовать OLS. .

person kmm    schedule 11.11.2011

arrow_upward
1
arrow_downward

В дополнение ко всем другим хорошим ответам вы можете использовать смещение. С категориальными предикторами немного сложнее, потому что вам нужно знать кодировку.

lm(weight~group+offset(1*(group=="Trt")))

1* здесь не нужен, но он вставлен, чтобы подчеркнуть, что вы проверяете гипотезу о том, что разница равна 1 (если вы хотите проверить гипотезу о разнице в d, используйте d*(group=="Trt")

person Ben Bolker    schedule 11.11.2011

arrow_upward
0
arrow_downward

Вы можете использовать t.test, чтобы сделать это для ваших данных. Параметр mu задает гипотезу о разнице групповых средних. Параметр alternative позволяет выбирать между односторонними и двусторонними тестами.

t.test(weight~group,var.equal=TRUE)

        Two Sample t-test

data:  weight by group 
t = 1.1913, df = 18, p-value = 0.249
alternative hypothesis: true difference in means is not equal to 0 
95 percent confidence interval:
 -0.2833003  1.0253003 
sample estimates:
mean in group Ctl mean in group Trt 
            5.032             4.661 



t.test(weight~group,var.equal=TRUE,mu=-1)

        Two Sample t-test

data:  weight by group 
t = 4.4022, df = 18, p-value = 0.0003438
alternative hypothesis: true difference in means is not equal to -1 
95 percent confidence interval:
 -0.2833003  1.0253003 
sample estimates:
mean in group Ctl mean in group Trt 
            5.032             4.661
person James    schedule 11.11.2011
comment
Хорошая альтернатива. Обратите внимание, что t-критерий доступен только для одной категориальной переменной, имеющей два уровня. - person kohske; 11.11.2011

arrow_upward
-1
arrow_downward

Закодируйте свой собственный тест. Вы знаете предполагаемый коэффициент и знаете стандартную ошибку. Вы можете создать свою собственную тестовую статистику.

person power    schedule 11.11.2011
comment
хотя ваш комментарий имеет смысл, вы должны помнить, что это R. поэтому вероятность того, что тривиальный тест уже реализован, чрезвычайно высока, и его стоит проверить, прежде чем тратить время на его кодирование. конечно, это имеет свою собственную обучающую ценность. - person Ramnath; 11.11.2011
comment
Это не обязательно плохой ответ, если вы добавили пример кода, чтобы проиллюстрировать, как это можно сделать. - person joran; 11.11.2011
comment
Это может быть тривиально для кода. Посмотрите учебник по эконометрике для первого года обучения. - person power; 11.11.2011