Sabatina,
assunto: Regularização
- Emerson Rigoni
- Henrique Aparecido Laureano [Lattes, LEG GitLab, GitLab]
Abril de 2016
Banco de dados
data("longley")Um dataset de 7 variáveis econômicas, observadas de 1947 até 1962 (n = 16)
GNP.deflator- deflator de preços implícito do GNP
GNP- Produto Nacional Bruto
Unemployed- Número de desempregados (1954 = 100)
Armed.Forces- número de pessoas nas forças armadas
Population- população ‘não institucionalizada’ \(\geq\) 14 anos de idade
Year- a ano (tempo)
Employed- número de pessoas empregadas
summary(longley) GNP.deflator GNP Unemployed Armed.Forces Population Year Employed
Min. : 83.00 Min. :234.3 Min. :187.0 Min. :145.6 Min. :107.6 Min. :1947 Min. :60.17
1st Qu.: 94.53 1st Qu.:317.9 1st Qu.:234.8 1st Qu.:229.8 1st Qu.:111.8 1st Qu.:1951 1st Qu.:62.71
Median :100.60 Median :381.4 Median :314.4 Median :271.8 Median :116.8 Median :1954 Median :65.50
Mean :101.68 Mean :387.7 Mean :319.3 Mean :260.7 Mean :117.4 Mean :1954 Mean :65.32
3rd Qu.:111.25 3rd Qu.:454.1 3rd Qu.:384.2 3rd Qu.:306.1 3rd Qu.:122.3 3rd Qu.:1958 3rd Qu.:68.29
Max. :116.90 Max. :554.9 Max. :480.6 Max. :359.4 Max. :130.1 Max. :1962 Max. :70.55 lattice:: splom(longley, pscales = 0, col = 1, type = c("p", "g", "smooth"), xlab = NULL
, main = "Todos os possíveis gráficos de dispersão 2 x 2"
, sub = "Curvas de tendência estimadas por suavização loess")
No gráfico acima podemos observar que várias variáveis são altamente correlacionadas, como por exemplo: Employed x Year, Year x Population, Employed x GNP, Year x GNP, Population x GNP, Year x GNP.deflator, Population x GNP.deflator, GNP x GNP.deflator
Regularização
library(glmnet)Por default a função glmnet utiliza um modelo linear Gaussiano
matriz de covariáveis
(x <- model.matrix(Employed ~ ., longley)[, -1]) GNP.deflator GNP Unemployed Armed.Forces Population Year
1947 83.0 234.289 235.6 159.0 107.608 1947
1948 88.5 259.426 232.5 145.6 108.632 1948
1949 88.2 258.054 368.2 161.6 109.773 1949
1950 89.5 284.599 335.1 165.0 110.929 1950
1951 96.2 328.975 209.9 309.9 112.075 1951
1952 98.1 346.999 193.2 359.4 113.270 1952
1953 99.0 365.385 187.0 354.7 115.094 1953
1954 100.0 363.112 357.8 335.0 116.219 1954
1955 101.2 397.469 290.4 304.8 117.388 1955
1956 104.6 419.180 282.2 285.7 118.734 1956
1957 108.4 442.769 293.6 279.8 120.445 1957
1958 110.8 444.546 468.1 263.7 121.950 1958
1959 112.6 482.704 381.3 255.2 123.366 1959
1960 114.2 502.601 393.1 251.4 125.368 1960
1961 115.7 518.173 480.6 257.2 127.852 1961
1962 116.9 554.894 400.7 282.7 130.081 1962vetor resposta
y <- longley$EmployedNa função glmnet temos o argumento alpha
alpha= 0, regressão ridgealpha= 1 (default), regressão lasso
Pra cada tipo de regressão (lasso e ridge) são ajustados três modelos, cada um com um \(\lambda\) diferente, e assim verificamos como isso impacta a estimação
Por default é usada uma sequência de 100 valores pra \(\lambda\), onde o maior valor é o menor valor para o qual todos os coeficientes são zero, e o menor é 0.0001, caso o número de observações seja maior que o número de variáveis, como acontece aqui
lasso
fit <- glmnet(x, y)
fit2 <- glmnet(x, y, lambda = seq(3, .001, length.out = 100))
fit3 <- glmnet(x, y, lambda = seq(3, .01, length.out = 100))
par(mfrow = c(1, 3))
plot(fit, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda default")
plot(fit2, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda: seq(3, 0.001, length.out = 100)")
plot(fit3, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda: seq(3, 0.01, length.out = 100)")
No gráfico da esquerda é usado o \(\lambda\) default
Na abscissa temos o log de \(\lambda\), na ordenada temos os valores do coeficientes
Cada curva é uma variável e no topo dos gráficos temos o número de coeficientes diferentes de zero
O maior log de \(\lambda\) é próximo de 1, logo, o maior \(\lambda\) é próximo de 3
3 é o maior \(\lambda\) considerado nos outros modelos
O número de \(\lambda\)s é fixo, 100
O que foi aumentando (leia-se, ficando mais próximo de zero) é o \(\lambda\) mínimo
No gráfico do meio ele é de 0.001 e no gráfico da direira é de 0.01, lembrando que o default é de 0.0001, ou seja, fomos aumentando o \(\lambda\) mínimo por uma razão de 10
Observações:
No \(\lambda\) default e no \(\lambda\) mínimo de 0.001 acabamos com todos os coeficientes diferentes de zero, quando aumentamos o \(\lambda\) mínimo pra 0.01 terminamos com 3 coeficientes de zero
ridge
fit4 <- glmnet(x, y, alpha = 0)
fit5 <- glmnet(x, y, alpha = 0, lambda = seq(3000, .001, length.out = 100))
fit6 <- glmnet(x, y, alpha = 0, lambda = seq(3000, .01, length.out = 100))
par(mfrow = c(1, 3))
plot(fit4, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda default")
plot(fit5, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda: seq(3000, 0.001, length.out = 100)")
plot(fit6, xvar = "lambda", ylab = "Coeficientes", las = 1, sub = "Lambda: seq(3000, 0.01, length.out = 100)")
O maior log de \(\lambda\) é próximo de 8, logo, o maior \(\lambda\) é próximo de 3000
3000 é o maior \(\lambda\) considerado nos outros modelos
Observações:
Para todos os \(\lambda\)s terminamos com todos os coeficientes diferentes de zero, na verdade os coeficientes se mostram diferentes de zero em todo o range da sequência de \(\lambda\)
E olhando para as estimativas dos coeficientes? (pra ambas regressões)
round(cbind(coefficients(fit, s = 0.0001), coefficients(fit4, s = 0.0001)), 5)7 x 2 sparse Matrix of class "dgCMatrix"
1 1
(Intercept) -2776.77868 -370.26268
GNP.deflator -0.02405 0.08373
GNP -0.00743 0.01072
Unemployed -0.01615 -0.00687
Armed.Forces -0.00932 -0.00166
Population -0.19185 0.12030
Year 1.47226 0.21049round(cbind(coefficients(fit2, s = 0.001), coefficients(fit5, s = 0.001)), 5)7 x 2 sparse Matrix of class "dgCMatrix"
1 1
(Intercept) -2546.62372 -2531.12692
GNP.deflator -0.00981 -0.01257
GNP -0.00638 -0.00550
Unemployed -0.01581 -0.01570
Armed.Forces -0.00899 -0.00897
Population -0.15003 -0.15417
Year 1.35094 1.34321round(cbind(coefficients(fit3, s = 0.01), coefficients(fit6, s = 0.01)), 5)7 x 2 sparse Matrix of class "dgCMatrix"
1 1
(Intercept) -1758.48727 -1345.39001
GNP.deflator . 0.03914
GNP . 0.00875
Unemployed -0.01377 -0.01274
Armed.Forces -0.00678 -0.00730
Population . -0.04652
Year 0.93628 0.72385Pegando o \(\lambda\) de 0.0001
Quando mudamos o tipo de regressão e deixamos o \(\lambda\) default as estimativas diferem consideravelmente! Maaaas alguns padrões ainda podem ser traçados, como por exemplo: O maior coeficiente ainda é o da variável Year, seguido de Population (com um detalhe, na regressão lasso o coeficiente é negativo, na regressão ridge ele é positivo!)
Com base nas estimativas podemos dizer que existe uma relação positiva com Year e uma relação negativa com Unemployed e Armed.Forces, i.e, conforme o ano passa o número de pessoas empregadas aumenta (Employed) e conforme o desemprego e o número de pessoas nas forças armadas aumenta, o número de pessoas desempregados diminui. Contudo, esses dois últimos coeficientes são pequenos, logo, essa relação não é muito grande
Pegando o \(\lambda\) de 0.001
Aqui as estimativas são muito mais próximas. O maior coeficiente ainda pertence ao Year e a menor a Popupation, i.e., conforme a população aumenta o número de pessoas empregadas diminui. Com exceção da variável Year, conforme as variáveis aumentam o valor da resposta diminui
Pegando o \(\lambda\) de 0.01
Na regressão lasso apenas três coeficientes foram diferentes de zero, Unemployed, Armed.Forces e Year. Na regressão ridge todos foram diferentes de zero. Nessas três variáveis existe uma certa conformidade de valores
Considerações:
A regressão ridge é conhecida por falhar na parcimônia do modelo, incluindo um número maior de variáveis no modelo. Isso pode justificar o comportamento observado com \(\lambda\) = 0.01
Em geral, quando não deixamos o \(\lambda\) no default os resultados obtidos com as duas regressões são similares, apontando para as mesmas direções
Lembrando que esse conjunto de dados é famoso por sua alta colinearidade
Com \(\lambda\)s maiores as regressão divergem em número de variáveis diferentes de zero, com \(\lambda\)s menores algumas estimativas se mostram diferentes e com \(\lambda\)s intermediários (\(\lambda\) de 0.001, nesse caso), as estimativas se mostram muito mais similares