Извод за коефициента на наклона в 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) и след това да погледнете само в посоката, определена от ast, но все още неизявената хипотеза. - 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