Exemplo Considere o modelo de regressão polinomial \[ Y = g(X) + \epsilon, \] onde \(X, Y\in \mathbb{R}\) e \(g(x) = \beta_0+\beta_1 x+\beta_2 x^2\).
Fornezida a amostra \((X_1,Y_1),\cdots,(X_n,Y_n)\) podemos estimar \(\beta= (\beta_0,\beta_1,\beta_2)\) com o estimador de mínimos quadrados \(\widehat{\beta}\). Suponha que \(g(x)\) seja côncava e estejamos interessados na localização do ponto \(x\) em que \(g(x)\) é maximizada.
Não é difícil encontrar que, o máximo ocorre em \(x =\theta\) onde \(\theta=-(1/2)\beta_1/\beta_2\). Uma estimativa pontual de \(\theta\) é \(\widehat{\theta} =-(1/2)\widehat{\beta}_1/\widehat{\beta}_2\).
Agora usamos o bootstrap para obter um intervalo de confiança para \(\theta\). A Figura 1 abaixo mostra 50 pontos extraídos do modelo acima com \(\beta_0 = 1\), \(\beta_1 = 2\), \(\beta_2 =-1\). Os \(X_i\)’s foram amostrados uniformemente em [0,2] e selecionamos \(\epsilon_i \sim N(0, 0.2^2)\), \(i=1,\cdots,50\).
Neste caso, \(\theta= 1\). As curvas verdadeira e estimada são mostradas na figura. Na parte inferior do gráfico, mostramos o intervalo de confiança do boostrap de 95% com base em \(B = 1.000\).
set.seed(74853)
beta.0 = -1; beta.1= 2; beta.2 = -1
eps = rnorm(50, mean = 0, sd = 0.2)
x = runif(50, min = 0, max = 2)
y = beta.0 + beta.1*x + beta.2*x^2 + eps
ajuste = lm(y ~ x+I(x^2))
summary(ajuste)
##
## Call:
## lm(formula = y ~ x + I(x^2))
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.3333 -0.1153 -0.0335 0.1070 0.5043
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.90992 0.09029 -10.078 2.49e-13 ***
## x 1.71628 0.21218 8.089 1.88e-10 ***
## I(x^2) -0.85866 0.10095 -8.506 4.53e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1749 on 47 degrees of freedom
## Multiple R-squared: 0.6096, Adjusted R-squared: 0.593
## F-statistic: 36.69 on 2 and 47 DF, p-value: 2.516e-10
theta = -(1/2)*coef(ajuste)[2]/coef(ajuste)[3]
theta
## x
## 0.9993968
library(bootstrap)
theta.func = function(x){
y = beta.0 + beta.1*x + beta.2*x^2 + eps
ajuste1 = lm(y ~ x+I(x^2))
return(-(1/2)*coef(ajuste1)[2]/coef(ajuste1)[3])
}
boot = bootstrap(x,1000,theta.func)
theta.ic <- function(x){mean(x)}
ic = boott(boot$thetastar,theta.ic)$confpoints
ic[3];ic[8]
## [1] 0.9995556
## [1] 1.002341
par(mar=c(3,4,1,1), pch = 19)
plot(x,y, xlim = c(0,2), ylim = c(-1,0.5), cex = 0.6)
xx = seq(0,2,by = 0.01)
lines(xx,beta.0+beta.1*xx+beta.2*xx^2, lwd = 2)
par(add=TRUE)
lines(xx,coef(ajuste)[1]+coef(ajuste)[2]*xx+coef(ajuste)[3]*xx^2, col = "red", lwd= 2)
abline(v=theta, lwd =2, col = "red")
grid()
Figura 1: 50 pontos retirados do modelo \(Y_i=-1+2X_i-X_i^2+\epsilon_i\), onde \(X_i\sim \mbox{Uniforme}[0.2]\) e \(\epsilon_i \sim N(0,0.2^2)\). Nesta situação o máximo do plinômio ocorre em \(\theta=1\). As curvas teórica e estimada são mostradas no gráfico.
Exercício No.1 Sejam \((X_1,Y_1,Z_1),\cdots,(X_n,Y_n,Z_n)\) vetores aleatórios independentes formados de triplas de números reais. O coeficientes de correlação parcial entre \(X\) e \(Y\) dado \(Z\) é definido como \[ \rho_{X,Y|Z} = -\frac{\Omega_{12}}{\sqrt{\Omega_{11}\Omega_{22}}}, \] onde \(\Omega=\Sigma^{-1}\) e \(\Sigma\) a matriz de covariancia de \((X,Y,Z)\).
O coeficiente de correlação parcial mede a dependência linear entre \(X\) e \(Y\) após remover o efeito de \(Z\). Para ilustração. suponha que geramos os dados da seguinte forma: escolhemos \(Z\sim N(0,1)\), \(X=10Z+\epsilon\) e \(Y=10Z+\delta\), onde \(\epsilon,\delta\sim N(0,1)\), independentes. A correlação entre \(X\) e \(Y\) é muito grande. Mas a correlação parcial é 0.
Encontrar o itervalo de confiança bootstrap com 95% de confiança para \(\rho_{X,Y|Z}\).