top of page

Regressão linear simples

Regressão linear múltipla

Diagnósticos

Regressão linear simples

​

Seja o seguinte banco de dados, que inclui a idade (age) e pressão arterial sistólica (SBP) de 30 indivíduos.

​

id age  SBP

 1  39  144

 2  47  220

 3  45  138

 4  47  145

 5  65  162

 6  46  142

 7  67  170

 8  42  124

 9  67  158

10  56  154

11  64  162

12  56  150

13  59  140

14  34  110

15  42  128

16  48  130

17  45  135

18  17  114

19  20  116

20  19  124

21  36  136

22  50  142

23  39  120

24  21  120

25  44  160

26  53  158

27  63  144

28  29  130

29  25  125

30  69  175

​

Fonte: http://people.sc.fsu.edu/~jburkardt/datasets/regression/x03.txt

 

dados<-read.table("clipboard",header=T)

​

# Ou então, use estas linhas para ler os dados:

​

urlfile <- "https://raw.githubusercontent.com/edsonzmartinez/cursoR/main/idadepressao.csv"

dados   <- read.csv2(urlfile,head=TRUE)

​

plot(dados$age,dados$SBP,xlab="Idade (anos)", ylab="Pressão arterial sistólica (mmHg)",bty="l",pch=19)

# Supondo que o valor de SBP para o indivíduo 2 é um erro de digitação

​

dados$age <- dados$age[dados$SBP<200]

​

dados$SBP <- dados$SBP[dados$SBP<200]

​

plot(dados$age,dados$SBP,xlab="Idade (anos) ",ylab="Pressão arterial sistólica (mmHg)",bty="l",pch=19)

# a função lm() (linear model) ajusta um modelo de regressão linear

​

lm(SBP~age,data=dados)

​

Call:

lm(formula = SBP ~ age)

 

Coefficients:

(Intercept)          age 

    97.0771       0.9493 

 

# a função summary() exibe mais detalhes sobre o ajuste do modelo

 

summary(lm(SBP~age,data=dados))

 

Call:

lm(formula = SBP ~ age)

 

Residuals:

    Min      1Q  Median      3Q     Max

-19.354  -4.797   1.254   4.747  21.153

 

Coefficients:

            Estimate Std. Error t value Pr(>|t|)   

(Intercept)  97.0771     5.5276  17.562 2.67e-16 ***

age           0.9493     0.1161   8.174 8.88e-09 ***

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

 

Residual standard error: 9.563 on 27 degrees of freedom

Multiple R-squared:  0.7122,    Adjusted R-squared:  0.7015

F-statistic: 66.81 on 1 and 27 DF,  p-value: 8.876e-09

 

# atribuindo ao objeto "modelo" os resultados do ajuste do modelo de regressão:

 

modelo <- lm(SBP~age,data=dados)

​

summary(modelo)

 

names(modelo)

​

 [1] "coefficients"  "residuals"     "effects"       "rank"        

 [5] "fitted.values" "assign"        "qr"            "df.residual" 

 [9] "xlevels"       "call"          "terms"         "model"

​

modelo$coefficients

​

(Intercept)         age

 97.0770843   0.9493225

 

class(modelo)

​

[1] "lm"

 

# Intervalos de confiança 95% para os parâmetros:

 

confint(modelo)

​

                 2.5 %     97.5 %

(Intercept) 85.7354850 108.418684

age          0.7110137   1.187631

 

# Gráfico:

​

plot(dados$age,dados$SBP,xlab="Idade (anos) ",ylab="Pressão arterial sistólica (mmHg)",bty="l",pch=19)

​

# Inserindo a reta de regressão no gráfico:

​

abline(modelo)

​

# Valores preditos:

 

fitted(modelo)

​

       1        2        3        4        5        6        7        8        9       10
134.1007 139.7966 141.6952 158.7830 140.7459 160.6817 136.9486 160.6817 150.2391 157.8337
     11       12       13       14       15       16       17       18       19       20
150.2391 153.0871 129.3541 136.9486 142.6446 139.7966 113.2156 116.0635 115.1142 131.2527
     21       22       23       24       25       26       27       28       29
144.5432 134.1007 117.0129 138.8473 147.3912 156.8844 124.6074 120.8101 162.5803

​

SBP


 [1] 144 138 145 162 142 170 124 158 154 162 150 140 110 128 130 135 114 116 124 136 142 120
[23] 120 160 158 144 130 125 175

 

# Gráfico de valores observados versus valores preditos:


plot(dados$SBP,fitted(modelo),ylab="Valores ajustados pelo modelo",xlab="Valores observados",pch=19,bty="l", main="Valores observados x preditos")

# Diagnósticos

 

par(mar = c(4, 4, 2, 2), mfrow = c(2, 2))

​

plot(modelo)

# Teste de normalidade dos resíduos:

​

shapiro.test(residuals(modelo))

 

        Shapiro-Wilk normality test

 

data:  residuals(modelo)

W = 0.96258, p-value = 0.38

 

# Tabela de análise de variância (ANOVA):

​

anova(modelo)

​

Analysis of Variance Table

 

Response: SBP

          Df Sum Sq Mean Sq F value    Pr(>F)   

age        1 6110.1  6110.1  66.808 8.876e-09 ***

Residuals 27 2469.3    91.5                     

---

Signif. codes: 

0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

 

# Gerando uma banda de confiança para a reta de regressão:

​

pconf <- predict(modelo,interval="confidence")

 

pconf

​

        fit      lwr      upr

1  134.1007 130.1804 138.0210

2  139.7966 136.1528 143.4404

3  141.6952 138.0225 145.3680

4  158.7830 152.7966 164.7695

5  140.7459 137.0954 144.3964

6  160.6817 154.3105 167.0529

7  136.9486 133.2322 140.6651

8  160.6817 154.3105 167.0529

...

​

datanew <- data.frame(cbind(dados$age,pconf))

names(datanew) <- c("age","fit","lwr","upr")

g <- datanew[with(datanew,order(dados$age)),]

plot(dados$age,dados$SBP,xlab="Idade (anos)",ylab="Pressão arterial sistólica (mmHg)",bty="l",pch=19)

points(g$age,g$lwr,type="l",col="red")

points(g$age,g$upr,type="l",col="red")

abline(modelo,col="blue")

#

# Alternativa:

 

plot(dados$age,dados$SBP,xlab="Idade (anos)",ylab="Pressão arterial sistólica (mmHg)",bty="l",pch=NA)

​

polygon(c(g$age,rev(g$age)),c(g$lwr,rev(g$upr)),col="gray",border=NA)

​

points(age,SBP,pch=19)

​

abline(modelo,col="blue")

​​

Regressão linear múltipla

 

Seja o banco de dados (clique aqui para abrir):

 

urlfile <- "https://raw.githubusercontent.com/edsonzmartinez/cursoR/main/altman.csv"

dados   <- read.csv2(urlfile,head=TRUE)

​

dados

​

   idade sexo altura  CTP

1     35    F    149 3.40

2     11    F    138 3.41

3     12    M    148 3.80

4     16    F    156 3.90

5     32    F    152 4.00

6     16    F    157 4.10

7     14    F    165 4.46

8     16    M    152 4.55

9     35    F    177 4.83

10    33    F    158 5.10

11    40    F    166 5.44

12    28    F    165 5.50

13    23    F    160 5.73

14    52    M    178 5.77

15    46    F    169 5.80

16    29    M    173 6.00

17    30    F    172 6.30

18    21    F    163 6.55

19    21    F    164 6.60

20    20    M    189 6.62

21    34    M    182 6.89

22    43    M    184 6.90

23    35    M    174 7.00

24    39    M    177 7.20

25    43    M    183 7.30

26    37    M    175 7.65

27    32    M    173 7.80

28    24    M    173 7.90

29    20    F    162 8.05

30    25    M    180 8.10

31    22    M    173 8.70

32    25    M    171 9.45

 

# A variável dependente é a CTP, capacidade total pulmonary (litros)

​

# A função levels() mostra os níveis de uma variável qualitativa 

​

levels(dados$sexo)
[1] "F" "M"

 

# A função contrasts() mostra como uma variável qualitativa é codificada na forma de variável dummy

​

contrasts(dados$sexo)
  M
F 0
M 1

​

# Regressão linear múltipla

 

modelo2 <- lm(CTP ~ idade+sexo+altura,data=dados)

 

modelo2

 

summary(modelo2)

 

Call:

lm(formula = CTP ~ idade + sexo + altura, data = dados)

 

Residuals:

    Min      1Q  Median      3Q     Max

-1.9572 -0.8404  0.0445  0.5793  2.6097

 

Coefficients:

            Estimate Std. Error t value Pr(>|t|)  

(Intercept) -8.54361    3.67861  -2.323  0.02770 *

idade       -0.02502    0.02353  -1.063  0.29665  

sexoM        0.69705    0.49944   1.396  0.17379  

altura       0.08955    0.02455   3.647  0.00107 **

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

 

Residual standard error: 1.156 on 28 degrees of freedom

Multiple R-squared:  0.5422,    Adjusted R-squared:  0.4932

F-statistic: 11.05 on 3 and 28 DF,  p-value: 5.822e-05

 

# Critério de Informação de Akaike (AIC)

​

AIC(modelo2)
[1] 105.808

​

# Diagnósticos

​

par(mar = c(4, 4, 2, 2), mfrow = c(2, 2))

​

plot(modelo2)

bottom of page