CE225 – Modelos Lineares Generalizadostaconeli/CE225/Aula1.pdf · 2015. 8. 5. · 2 Objetivo da...

1

CE225 – Modelos Lineares Generalizados

2

Objetivo da disciplina

• Apresentar ao aluno a teoria e aplicações dos Modelos Lineares Generalizados,

propostos originalmente em Nelder e Wedderburn (1972), que configuram extensões dos

modelos lineares clássicos (com erros normalmente distribuídos) e que permitem

analisar a relação funcional entre um conjunto de variáveis independentes e uma variável

aleatória dependente com distribuição pertencente à família exponencial de

distribuições.

• A família exponencial contempla, dentre outras, as distribuições normal, exponencial,

gama, normal inversa, Poisson, binomial e binomial negativa.

• Na sequência são descritos alguns dos problemas que serão analisados, ao longo do

semestre, usando elementos de Modelos Lineares Generalizados.

3

Exemplos de motivação

4

Exemplo 1 – Análise da resistência de uma nova fibra sintética usada na produção de

camisas. Sabendo-se que a resistência da fibra é afetada pela quantidade de algodão

utilizada, e que a quantidade de algodão no produto final, de acordo com as características

desejadas, deve estar no intervalo de 10 a 40%, um experimento é delineado com cinco

réplicas (amostras de tecidos) para cinco diferentes especificações referentes à quantidade de

algodão.

Variável resposta: Resistência da fibra (em libras/pol2).

Variável explicativa: Porcentagem de algodão no tecido, com cinco níveis: 15, 20, 25, 30 e

35%.

5

Dados

Quadro 1 – Dados de resistência (em libras/pol2) para o experimento

de fibra sintética.

Amostra de tecido Porcentagem

de algodão 1 2 3 4 5

15 7 7 15 11 9

20 12 17 12 18 18

25 14 18 18 19 19

30 19 25 22 19 23

35 7 10 11 15 11

Objetivos:

• Avaliar o efeito da porcentagem de algodão na resistência da fibra sintética;

• Identificar a porcentagem ideal de algodão de forma a se obter máxima resistência.

6

15 20 25 30 35

10

15

20

25

Porcentagem de algodão

Res

istê

ncia

da

fibra

Figura 1 – Gráfico de dispersão para as resistências das fibras sob cinco porcentagens distintas de

algodão.

7

Exemplo 2 – Amostras de 20 insetos (Heliothis virescens - praga do algodão) foram expostas

a doses crescentes do cipermetrina, dois dias depois da emergência da pupa (Collet, 2002).

Após 72h, foram contados os insetos mortos.

Variável resposta: Número de insetos mortos.

Variáveis explicativas

• Dose de cipermetrina, com níveis 1, 2, 4, 8, 16, 32 u.m.:

• Sexo (macho ou fêmea).

8

Quadro 2 – Números de insetos mortos

em amostras de 20 insetos machos e fêmeas

submetidos a doses crescentes de

cipermetrina.

Nº insetos mortos Dose Log(Dose)

Machos Fêmeas

1,0 0 1 0

2,0 1 4 2

4,0 2 9 6

8,0 3 13 10

16,0 4 18 12

32,0 5 20 16

9

Objetivos:

• Propor um modelo que descreva o aumento na mortalidade dos insetos segundo a dose de

aplicada de cipermetrina;

• Comparar as curvas de mortalidade para insetos machos e fêmeas;

• Estimar doses letais, ou seja, doses efetivas (mortais) para uma proporção p de insetos.

10

0 1 2 3 4 5

0.0

0.2

0.4

0.6

0.8

1.0

log(dose)

Pro

porç

ão d

e in

seto

s m

orto

s

Machos

Fêmeas

Figura 2 – Gráfico da proporção de insetos mortos segundo o sexo e a dose de inseticida.

11

Exemplo 3 - Uma população de mulheres indígenas vivendo numa região próxima a Phoenix,

Arizona, foi testada para diabetes de acordo com o critério estabelecido pela Organização Mundial

de Saúde. Os dados foram coletados pelo Instituto Nacional de Diabetes e Doenças Digestivas e

Renais dos EUA. São considerados os dados referentes aos 532 registros completos.

Variáveis explicativas

npreg - número de gestações;

gli - concentração de glicose no plasma no teste de tolerância à glicose oral.

bp - pressão arterial diastólica (mm Hg).

skin - espessura da prega tricipital (mm).

bmi - índice de massa corporal (peso/altura2).

ped - função pedigree diabetes.

age - idade em anos.

Variável resposta

type – diagnóstico de diabetes de acordo com o teste de Glicemia em Jejum.

12

Notas:

O Teste Oral de Tolerância à Glicose (também conhecido como Curva Glicêmica) é feito da

seguinte maneira: a pessoa com suspeita de diabetes ingere 75g de glicose diluída em água. Após duas

horas de espera, é feita a coleta de sangue para medir a taxa de glicose. Se o resultado for igual ou

superior a 200mg/dl (miligramas por decilitro), o indivíduo é considerado portador de diabetes. Se a

glicemia estiver entre 140mg/dl e 199mg/dl, então o diagnóstico é de pré-diabetes.

13

Dados - Seis primeiras linhas da base:

npreg glu bp skin bmi ped age type

1 5 86 68 28 30.2 0.364 24 No

2 7 195 70 33 25.1 0.163 55 Yes

3 5 77 82 41 35.8 0.156 35 No

4 0 165 76 43 47.9 0.259 26 No

5 0 107 60 25 26.4 0.133 23 No

6 5 97 76 27 35.6 0.378 52 Yes

Objetivos:

• Identificar fatores associados à incidência de diabetes;

• Estabelecer um modelo preditivo para o diagnóstico de diabetes.

14

Não Sim

0

2

4

6

8

10

12

14

Diagnóstico

Gra

vide

z

Não Sim

60

80

100

120

140

160

180

200

DiagnósticoG

licos

e

Não Sim

40

50

60

70

80

90

100

110

Diagnóstico

Pre

ssão

art

eria

l (m

mH

g)

Não Sim

20

40

60

80

100

Diagnóstico

Esp

essu

ra p

rega

tri

c.

Não Sim

20

25

30

35

40

45

Diagnóstico

IMC

(kg

/m2)

Não Sim

0.0

0.5

1.0

1.5

2.0

Diagnóstico

Ped

igre

e

Não Sim

20

30

40

50

60

Diagnóstico

Idad

e (a

mos

)

Não Sim

Diagnóstico

Núm

ero

de p

acie

ntes

0

20

40

60

80

100

120

Figura 3 – Distribuição das variáveis explicativas segundo o diagnóstico de diabetes.

15

Tabela 1 – Médias, desvios padrões e estatística do teste t

(comparação de duas médias), para amostras independentes,

para as variáveis explicativas segundo o diagnóstico

Diagnóstico Variável

Não Sim

Estatística

t

Número de gestações 2,9 (2,8) 4,8 (4,0) -3,56

Glicose (oral) 113 (26) 145 (30) -7,38

Pressão diastólica 69 (11) 74 (11) -2,95

Espessura – prega tricipital 27 (11) 33 (12) -3,39

Pedigree 0,41 (0,27) 0,55 (0,36) -4,51

Idade 29 (10) 37 (11) -2,70

16

npreg

60 100 140 180 20 40 60 80 100 0.0 0.5 1.0 1.5 2.0

02468101214

6080

100120140160180200

glu

bp405060708090100110

20406080

100

skin

bmi202530354045

0.00.51.01.52.0

ped

0 2 4 6 8 12 40 60 80 100 20 30 40 20 30 40 50 60

20

3040

5060

age

Figura 4 – Matriz de gráficos de dispersão para as variáveis explicativas.

17

Exemplo 4 – Na sequência, são apresentadas as cinco primeiras linhas de um banco de dados

referente a um estudo prospectivo com 100 indivíduos de pelo menos 65 anos de idade em boas

condições físicas. O objetivo do estudo é tentar relacionar o número de quedas num período de

seis meses com as seguintes variáveis explicativas, descritas na ordem em que aparecem na base:

Variáveis explicativas:

• Intervenção – Fator com níveis ‘Educ’: educação somente; ‘Educ+Exerc’: educação e exercícios

físicos;

• Sexo – Fator com níveis ‘Fem’: feminino; ‘Masc’: masculino;

• Balanço – escore de equilíbrio do indivíduo, numa escala de 0 a 100 (quanto maior o escore,

maior o equilíbrio;

• Força – escore de força do indivíduo, numa escala de 0 a 100 (quanto maior o escore, maior a

força).

Variável resposta:

• Quedas – número de quedas no período;

18

Objetivos:

• Avaliar o efeito da intervenção na prevenção das quedas;

• Identificar características dos indivíduos associadas a um maior número de quedas.

19

Dados: 10 primeiras linhas da base:

> head(geriatra,10)

quedas intervenção sexo balanço força

1 1 Educ+Exerc Fem 45 70

2 1 Educ+Exerc Fem 62 66

3 2 Educ+Exerc Masc 43 64

4 0 Educ+Exerc Masc 76 48

5 2 Educ+Exerc Fem 51 72

6 1 Educ+Exerc Masc 73 39

7 0 Educ+Exerc Masc 40 54

8 0 Educ+Exerc Fem 66 37

9 2 Educ+Exerc Masc 80 81

10 2 Educ+Exerc Masc 56 60

20

• Análise descritiva (univariada)

> summary(geriatra)

quedas intervenção sexo balan ço força

Min. : 0.00 Educ :50 Fem :47 Min. : 13.00 Min. :18.00

1st Qu.: 1.00 Educ+Exerc:50 Masc:53 1st Qu.: 39.00 1st Qu.:52.00

Median : 3.00 Median : 51.50 Median :60.00

Mean : 3.04 Mean : 52.83 Mean :60.78

3rd Qu.: 4.00 3rd Qu.: 66.25 3rd Qu.:70.25

Max. :11.00 Max. : 98.00 Max. :90.00

21

• Análise descritiva (bivariada):

o Número de quedas vs intervenção;

> with(geriatra,describeBy(quedas, intervenção, mat = TRUE,digits=2))

item group1 vars n mean sd median trimmed mad min max range skew kurtosis se

11 1 Educ 1 50 4.52 2.40 4 4.25 1.48 1 11 10 0.89 -0.07 0.34

12 2 Educ+Exerc 1 50 1.56 1.33 1 1.43 1.48 0 5 5 0.62 -0.52 0.19

o Número de quedas vs sexo;

> with(geriatra,describeBy(quedas, sexo, mat = TRUE ,digits=2))

item group1 vars n mean sd median trimmed ma d min max range skew kurtosis se

11 1 Fem 1 47 3.47 2.49 3 3.21 1.4 8 0 11 11 0.98 0.70 0.36

12 2 Masc 1 53 2.66 2.34 2 2.35 1.4 8 0 10 10 1.18 1.08 0.32

22

o Número de quedas vs nível de balanço;

> with(geriatra,describeBy(quedas, cut(balanço,4), mat = TRUE,digits=2))

item group1 vars n mean sd median trimme d mad min max range skew kurtosis se

11 1 (12.9,34.2] 1 19 2.32 1.25 2 2.3 5 1.48 0 4 4 -0.09 -1.34 0.29

12 2 (34.2,55.5] 1 39 2.95 2.70 2 2.5 8 1.48 0 11 11 1.36 1.35 0.43

13 3 (55.5,76.8] 1 31 3.32 2.60 3 3.1 6 2.97 0 9 9 0.43 -0.92 0.47

14 4 (76.8,98.1] 1 11 3.82 2.44 3 3.6 7 1.48 0 9 9 0.68 -0.31 0.74

o Número de quedas vs nível de força;

> with(geriatra,describeBy(quedas, cut(força,4), ma t = TRUE,digits=2))

item group1 vars n mean sd median trimmed mad min max range skew kurtosis se

11 1 (17.9,36] 1 5 2.80 1.30 3 2.80 1.48 1 4 3 -0.26 -1.96 0.58

12 2 (36,54] 1 25 2.32 2.19 2 2.10 2.97 0 7 7 0.61 -0.74 0.44

13 3 (54,72] 1 50 3.30 2.55 3 2.92 1.48 0 10 10 1.10 0.25 0.36

14 4 (72,90.1] 1 20 3.35 2.60 3 3.06 1.48 0 11 11 1.10 1.45 0.58

23

0 1 2 3 4 5

-4-2

02

46

Fitted values

Res

idua

ls

Residuals vs Fitted

5293

67

-2 -1 0 1 2

-10

12

3

Theoretical Quantiles

Sta

ndar

dize

d re

sidu

als

Normal Q-Q

52

93

67

0 1 2 3 4 5

0.0

0.5

1.0

1.5

Fitted values

Sta

ndar

dize

d re

sidu

als

Scale-Location52

93

67

0.00 0.05 0.10 0.15-2

-10

12

3

Leverage

Sta

ndar

dize

d re

sidu

als

Cook's distance

0.5

Residuals vs Leverage

5293

42

Figura 5 – Gráficos para os resíduos de um modelo de regressão linear múltipla ajustado aos dados de

quedas de idosos.

24

Exemplo 5 – Dados de 4624 apólices de seguros de automóveis que registraram sinistros no

período de um ano, entre 2004 e 2005.

Variáveis explicativas:

• Valor: Valor do veículo (x10.000 dólares);

• Tipo: tipo de veículo, com níveis: o BUS: ônibus; o CONVT: conversível; o COUPE; o HBACK: hathback; o HDTOP: hardtop; o MCARA: trailer motorizado; o MIBUS: mini-ônibus; o RDSTR: roadster; o SEDAN; o STNWG: station wagon; o TRUCK: caminhão; o UTE: utilitário.

• Idade: idade do veículo, com níveis 1 (veículos mais novos), 2, 3 e 4;

• Sexo: sexo do motorista, com níveis: M (masculino) e F (Feminino);

• Area: área de residência do motorista: A, B, C, D, E e F;

• Idademot: idade do motorista, com níveis: 1 (mais novos), 2, 3, 4, 5 e 6. Variável resposta:

• Quantia: valor (somado) dos sinistros apresentados no período (variável resposta).

25

Objetivos:

• Identificar fatores associados a sinistros mais caros;

• Estabelecer um modelo preditivo, que sirva de base para estabelecer a tabela de preços para

períodos futuros.

26

Dados

• Dez primeiras linhas da base:

> head(dados2,10)

Valorcar Quantia Tipo Idade Sexo Area Idadem ot

15 1.66 669.5100 SEDAN 3 M B 6

17 1.51 806.6100 SEDAN 3 F F 4

18 0.76 401.8055 HBACK 3 M C 4

41 1.89 1811.7100 STNWG 3 M F 2

65 4.06 5434.4400 STNWG 2 M F 3

66 1.39 865.7900 HBACK 3 F A 4

96 2.66 1105.7700 STNWG 1 F F 5

99 0.50 200.0000 HBACK 4 F A 5

116 1.16 739.2300 STNWG 4 F B 2

125 3.56 3230.6000 MCARA 3 M F 4

27

• Análise descritiva (univariada):

> summary(dados2)

Valorcar Quantia Tipo Idade Sexo Area Idademot

Min. : 0.000 Min. : 200.0 SEDAN :1476 1: 825 F:2648 A:1085 1: 496

1st Qu.: 1.100 1st Qu.: 353.8 HBACK :1264 2:1259 M:1976 B: 965 2: 932

Median : 1.570 Median : 761.6 STNWG :1173 3:1362 C:1412 3:1113

Mean : 1.859 Mean : 2014.4 UTE : 260 4:1178 D: 496 4:1104

3rd Qu.: 2.310 3rd Qu.: 2091.4 HDTOP : 130 E: 386 5: 614

Max. :13.900 Max. :55922.1 TRUCK : 120 F: 280 6: 365

(Other): 201

> sum(dados2$Quantia>15000) ### Numero de apólices que geraram mais de

$15.000 em sinistros.

[1] 65

28

Todas as apólices com sinistros

Valor dos sinistros (x$10.000)

Fre

quên

cia

0 5 10 15 20 25 30 35 40 45 50 55 60

0

200

400

600

800

1000

1200

1400

1600

1800

Apólices com sinistros inferiores a $15.000

Valor dos sinistros (x$10.000)

Fre

quên

cia

0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

0

200

400

600

800

1000

1200

1400

1600

1800

Figura 6 – Distribuição dos valores dos sinistros gerados pelos segurados.

29

• Análise descritiva (bivariada) – Valores dos sinistros vs tipo de veículo:

> with(dados2,describeBy(Quantia, Tipo, mat = TRUE,digits=2))

item group1 vars n mean sd median trimmed mad min max range skew kurtosis se

11 1 BUS 1 9 1484.79 1483.53 876.48 1484.79 613.99 371.82 4790.84 4419.02 1.13 -0.11 494.51

12 2 CONVT 1 3 2296.27 3319.80 530.00 2296.27 440.33 233.00 6125.81 5892.81 0.38 -2.33 1916.69

13 3 COUPE 1 68 2760.64 4197.29 1171.92 1819.55 1440.97 200.00 19847.74 19647.74 2.66 7.15 509.00

14 4 HBACK 1 1264 2048.37 3291.24 783.58 1313.14 865.22 200.00 47296.61 47096.61 4.52 37.81 92.57

15 5 HDTOP 1 130 2267.78 5015.92 584.16 1064.98 484.53 200.00 32814.80 32614.80 4.26 19.88 439.93

16 6 MCARA 1 14 762.42 812.92 379.97 591.53 38.81 345.00 3230.60 2885.60 2.02 3.21 217.26

17 7 MIBUS 1 43 2700.11 4529.85 1286.59 1471.73 1068.48 200.00 20545.10 20345.10 2.67 6.29 690.80

18 8 PANVN 1 62 2146.99 3552.48 714.79 1333.94 626.92 200.00 22216.09 22016.09 3.41 14.50 451.17

19 9 RDSTR 1 2 684.73 685.51 684.73 684.73 718.66 200.00 1169.46 969.46 0.00 -2.75 484.73

110 10 SEDAN 1 1476 1816.82 2928.87 759.81 1144.24 829.98 200.00 29634.63 29434.63 4.09 24.18 76.24

111 11 STNWG 1 1173 2014.57 4063.93 734.19 1164.22 674.26 200.00 55922.13 55722.13 6.27 55.03 118.66

112 12 TRUCK 1 120 2662.47 4675.85 807.55 1392.42 672.78 200.00 22405.44 22205.44 2.81 7.30 426.84

113 13 UTE 1 260 2296.96 3728.77 782.51 1384.33 755.53 200.00 28012.83 27812.83 3.35 14.33 231.25

30

• Valores dos sinistros vs idade do motorista:

> with(dados2,describeBy(Quantia, Idademot, mat = TRUE,digits=2))

item group1 vars n mean sd median trimmed mad min max range skew kurtosis se

11 1 1 1 496 2635.83 4320.54 994.66 1705.46 1178.16 200 46868.18 46668.18 4.37 29.93 194.00

12 2 2 1 932 2129.66 4106.12 760.19 1241.36 830.55 200 55922.13 55722.13 5.64 47.35 134.50

13 3 3 1 1113 1915.64 3065.13 743.77 1201.50 806.19 200 31974.77 31774.77 3.74 19.19 91.88

14 4 4 1 1104 1943.21 3503.93 750.47 1165.36 816.13 200 47296.61 47096.61 5.36 43.88 105.46

15 5 5 1 614 1728.68 2798.37 702.76 1054.16 745.40 200 22216.09 22016.09 3.68 16.95 112.93

16 6 6 1 365 1872.79 3405.27 852.32 1110.64 826.48 200 31243.67 31043.67 4.69 27.66 178.24

31

BUS

CONVT

COUPE

HBACK

HDTOP

MCARA

MIBUS

PANVN

RDSTR

SEDAN

STNWG

TRUCK

UTE

0

10

20

30

40

50

Valor dos sinistros x $10.000

1

2

3

4

0

10

20

30

40

50

Idade d

o m

oto

rista

Valor dos sinistrosx $10.000

F

M

0

10

20

30

40

50

Sexo

Valor dos sinistrosx $10.000

A

B

C

D

E

F

0

10

20

30

40

50

Áre

a d

e re

sid

ência

Valor dos sinistrosx $10.000

1

2

3

4

5

6

0

10

20

30

40

50

Idade d

o m

oto

rista

Valor dos sinistrosx $10.000

0

2

4

6

8

10

12

14

0

10

20

30

40

50

Valo

r do ve

ículo

Valor dos sinistrosx $10.000

Fig

ur

a 7

– D

istribu

ição

do

s va

lores d

os sin

istros seg

un

do

as co

va

riáv

eis (tod

as a

s ap

ólices co

m sin

istro).

32

Exemplo 6 – Análise do desempenho de cinco tipos de turbinas de aviões. Foi conduzido um

experimento com 10 turbinas de cada tipo

Variável explicativa – Tipo de turbina, com níveis I, II, III, IV e V.

Variável resposta – Tempo de vida da turbina, em milhões de ciclos até verificada a perda de

velocidade.

Objetivos –

• Estimar parâmetros correspondentes às distribuições dos tempos de vida dos cinco tipos de

turbinas;

• Comparar os tempos médios de vida, identificar quais turbinas são mais resistentes.

33

Dados:

> turbdata

I II III IV V

1 3.03 3.19 3.46 5.88 6.43

2 5.53 4.26 5.22 6.74 9.97

3 5.60 4.47 5.69 6.90 10.39

4 9.30 4.53 6.54 6.98 13.55

5 9.92 4.67 9.16 7.21 14.45

6 12.51 4.69 9.40 8.14 14.72

7 12.95 5.78 10.19 8.59 16.81

8 15.21 6.79 10.71 9.80 18.39

9 16.04 9.37 12.58 12.28 20.84

10 16.84 12.75 13.41 25.46 21.51

34

Análise descritiva:

Turbina 1 Turbina 2 Turbina 3 Turbina 4 Turbina 5

510

1520

25

Tipo de turbina

Tem

po a

té p

erda

de

velo

cida

de (m

ilhõe

s de

cic

los)

Figura – Distribuição dos tempos de vida segundo tipo de turbina.

35

> medias=with(turbdata,tapply(Tempo,Turbina,mean)); medias

Turbina 1 Turbina 2 Turbina 3 Turbina 4 Turbina 5

10.693 6.050 8.636 9.798 14.706

> variancias=with(turbdata,tapply(Tempo,Turbina,var )); variancias

Turbina 1 Turbina 2 Turbina 3 Turbina 4 Turbina 5

23.225512 8.497489 10.828116 33.711796 23.652316

> cvs=sqrt(variancias)/medias;cvs # Coeficientes de variação.

Turbina 1 Turbina 2 Turbina 3 Turbina 4 Turbina 5

0.4506954 0.4818257 0.3810341 0.5925889 0.3307061