R: Práctica de Análisis de Regresión con Dos Variables

Análisis de Regresión con Dos Variables

A continuación encontrarán el script que se ha desarrollado durante las clases de Econometría 1 para el tema de Análisis de Regresión con Dos Variables. Puede copiar y pegar este script en un editor de R para correr las instrucciones junto a las notas de clases.

#----Econometría 1 - Prof. Anthonny Arias----#

#--Limpiamos nuestro espacio de trabajo--#

rm(); rm(list=ls())
cat("\014")

# Definimos la variable escolaridad y su media.
# Para esto, usamos la instrucción c() para definir vector.
# Y usamos la instrucción mean() para definir la media.

escolaridad <- c(6,7,8,9,10,11,12,13,14,15,16,17,18)
m.escolaridad <- mean(escolaridad)

# Definimos la variable salario y su media.

salario <- c(4.4567,5.77,5.9787,7.3317,7.3182,6.5844,7.8182,7.8351,11.0223,10.6738,10.8361,13.615,13.531)
m.salario <- mean(salario)

# Hacemos un gráfico de dispersión de estas dos variables.

plot(escolaridad,salario)

# Una vez obtenidos estos valores, podemos calcular los estimadores beta1 y beta2.

beta2 <- sum( (escolaridad-m.escolaridad)*(salario-m.salario) )/sum( (escolaridad-m.escolaridad)^2 )
beta2

beta1 <- m.salario - beta2*m.escolaridad
beta1

# Calculamos los valores estimados del salario.

salario.e <- beta1 + beta2*escolaridad
salario.e

# Calculamos los residuos.

residuos <- salario - salario.e
residuos

# Calculamos la var.e.

var.e <- sum( (residuos)^2 )/(length(salario)-2)

var.e
# Caculamos el error estándar, aplicando la raíz cuadrada a la var.e.

error.s <- sqrt(var.e)
error.s

# Calculamos la var.e de beta2

v.beta2 <- var.e/sum( (escolaridad-m.escolaridad)^2 )
v.beta2

# Calculamos el error estándar de beta2

es.beta2 <- sqrt(v.beta2)
es.beta2

# Para calcular el intervalo de confianza de beta2, consideramos t=1.7959

li.beta2 <- beta2 - qt(0.975,df=length(escolaridad)-2)*es.beta2
li.beta2

ls.beta2 <- beta2 + qt(0.975,df=length(escolaridad)-2)*es.beta2
ls.beta2

# Calculamos la var.e de beta1.

v.beta1 <- var.e*sum( escolaridad^2 )/(length(escolaridad) * sum( (escolaridad-m.escolaridad)^2 ))
v.beta1

# Calculamos el error estándar de beta1

es.beta1 <- sqrt(v.beta1)
es.beta1

# Para calcular el intervalo de confianza de beta1, consideramos t=1.7959

li.beta1 <- beta1 - qt(0.975,df=length(escolaridad)-2)*es.beta1
li.beta1

ls.beta1 <- beta1 + qt(0.975,df=length(escolaridad)-2)*es.beta1
ls.beta1

# Para hacer la prueba de hipótesis bilateral, determinamos el t-calculado.
t.c <- (beta2-0.70)/es.beta2
t.c # Como t.c está fuera del intervalo (-2.201,2.201) entonces rechazamos la hipótesis nula.
qt(0.025,df=length(escolaridad)-2)
qt(0.975,df=length(escolaridad)-2)


# Para hacer la prueba de hipótesis unilateral, determinamos el t-calculado.
t.c <- (beta2-0.50)/es.beta2
t.c # Como t.c está fuera del intervalo (-2.201,2.201) entonces rechazamos la hipótesis nula.
qt(0.95,df=length(escolaridad)-2)

# Calculamos ahora, el intervalo de confianza para chi-cuadrado
li.var.e <- (length(escolaridad)-2)*var.e/qchisq(0.975,df=11)
li.var.e

ls.var.e <- (length(escolaridad)-2)*var.e/qchisq(0.025,df=11)
ls.var.e

# Como la hipótesis nula indica que la varianza es igual a 0.6, entonces no rechazamos esta hipótesis.

# Podemos también llevar a cabo esta prueba con el estadístico chi-cuadrado. Para esto, calculamos el estadístico chi-cuadrado.

chi.c <- (length(escolaridad)-2)*var.e/0.6
chi.c

li.chi <- qchisq(0.025,df=df=length(escolaridad)-2)
li.chi
ls.chi <- qchisq(0.975,df=df=length(escolaridad)-2)
ls.chi

# Éste está dentro del intervalo [ qchisq(0.025,df=df=length(escolaridad)-2) ; qchisq(0.025,df=df=length(escolaridad)-2) ], por lo tanto, no se rechaza H0.

# Calculamos la suma de los cuadrados explicada.

SCE.escolaridad <- beta2^2*sum( (escolaridad-m.escolaridad)^2 )
SCP.escolaridad <- SCE.escolaridad/1
SCP.escolaridad

# Calculamos la suma de los cuadrados de los residuos.

SCR.residuos <- sum(residuos^2)
SCP.residuos <- SCR.residuos/(length(escolaridad)-2)
SCP.residuos

# Calculamos la suma de los cuadrados totales.

SCT.salarios <- sum( (salario-m.salario)^2 )
SCP.salarios <- SCT.salarios/(length(salario)-1)
SCP.salarios

# Calculamos ahora el valor F (F calculado).

F.c <- SCP.escolaridad/SCP.residuos
F.c

# Calculamos el p-value (valor-p) para este F calculado.

pf(F.c,1,length(escolaridad)-2,lower.tail = F)

# Verificamos que se cumple el teorema

t.c <- (beta2-0)/es.beta2
t.c
t.c^2
F.c

# Predicción de la Media

escolaridad.0 <- 20
salario.0 <- beta1 + beta2*escolaridad.0
salario.0

# Calculamos la varianza de la predicción.

varm.salario.0 <- var.e*(1/length(escolaridad)+(escolaridad.0-m.escolaridad)^2/sum((escolaridad-m.escolaridad)^2))
varm.salario.0

# Calculamos ahora el error estándar.

eem.salario.0<- sqrt(varm.salario.0)
eem.salario.0

# Calculamos el intervalo de confianza para salario.0

li.salario.0 <- beta1 + beta2*escolaridad.0 - qt(0.025,df=length(escolaridad)-2,lower.tail = FALSE)*eem.salario.0
li.salario.0

ls.salario.0 <- beta1 + beta2*escolaridad.0 + qt(0.025,df=length(escolaridad)-2,lower.tail = FALSE)*eem.salario.0
ls.salario.0

# Predicción Individual

# Calculamos la varianza de la predicción.

vari.salario.0 <- var.e*(1+1/length(escolaridad)+(escolaridad.0-m.escolaridad)^2/sum((escolaridad-m.escolaridad)^2))
vari.salario.0

# Calculamos ahora el error estándar.

eei.salario.0<- sqrt(vari.salario.0)
eei.salario.0

# Calculamos el intervalo de confianza para salario.0

li.salario.0 <- beta1 + beta2*escolaridad.0 - qt(0.025,df=length(escolaridad)-2,lower.tail = FALSE)*eei.salario.0
li.salario.0

ls.salario.0 <- beta1 + beta2*escolaridad.0 + qt(0.025,df=length(escolaridad)-2,lower.tail = FALSE)*eei.salario.0
ls.salario.0

#----Análisis de Residuos----#

#--Análisis de Correlación--#

# Gráfico de dispersión para los residuos.

plot(residuos)

# Hacemos la gráfica de la función de autocorrelación.
# Si todaslas barras están por debajo de las líneas azules, esto indica que no hay autocorrelación.
# https://www.reddit.com/r/AskStatistics/comments/5kiix2/interpret_acfpacf_dataplots_in_r/

acf(residuos)

# Hacemos la prueba de Durbin–Watson, que establece como hipótesis nula que el coeficiente de correlación es igual a cero.
# El estadístico de Durbin–Watson igual a dos indica que no hay autocorrelación.
# https://en.wikipedia.org/wiki/Durbin%E2%80%93Watson_statistic

library("lmtest")
dwtest(salario ~ escolaridad)

#--Pruebas de Normalidad--#

# Generamos el histograma de los residuos.

hist(residuos)

plot(density(residuos))

# Gráfica de probabilidad normal

qqnorm(residuos, pch = 1, frame = TRUE)
qqline(residuos, col = "steelblue", lwd = 2)

# También se puede llevar a cabo usando el siguiente comando
library("car")
qqPlot(residuos,col.lines="steelblue")

# Prueba de Anderson-Darling.

library(nortest)
ad.test(residuos)

# Prueba de normalidad de Jarque-Bera (JB)
# https://lancebachmeier.com/computing/j-b-test.html
# Esta plantea como hipótesis nula el coeficiente de asimetría igual cero y la curtosis igual a tres.

library(tseries)
jarque.bera.test(residuos)

# Prueba de normalidad de Shapiro-Wilk
# https://stat.ethz.ch/R-manual/R-devel/library/stats/html/shapiro.test.html
# Esta prueba plantea como hipótesis nula que los datos están normalmente distribuídos.

shapiro.test(residuos)

R: Estimadores Mínimos Cuadrados Ordinarios (MCO)

El análisis de regresión sienta la base para los estudios econométricos y a su vez, estos se fundamentan formulando modelos lineales con dos variables: una independiente y otra dependiente; este tipo de modelos definen rectas, es decir, aquellos que se expresan de la siguiente forma:

Y = \beta_1 + \beta_2 X

También pudiera interesarte

Anuncios

Linealidad

Al mencionar la linealidad en una relación entre variables, siempre es importante especificar respecto a qué elemento de la relación, es dicha relación, lineal. Formalmente, diremos que una relación es lineal respecto a un elemento de la ecuación, si dicho elemento no está siendo multiplicada por sí mismo o si permanece inalterado por alguna función en la expresión, por ejemplo, la siguiente ecuación

Y = \beta_1^2 + \beta_2 \cdot \ln(X)

Es una ecuación lineal respecto respecto la variable Y y el parámetro \beta_2, debido a que estos dos elementos permanecen inalterados. Sin embargo, no es lineal respecto al parámetro \beta_1 pues este está multiplicado por sí mismo, tampoco es lineal respecto a la variable X pues esta está alterada por la función logaritmo neperiano.

La linealidad respecto a los parámetros representa una base en la que se fundamentan los Modelos Lineales que estudiaremos. Es por esto que, usualmente, el término regresión lineal hace referencia a la linealidad de los parámetros. Por lo tanto, puede o no ser lineal en las variables.

El Modelo de Regresión Lineal

Todo estudio de índole estadístico está sometido a un error de aproximación y la econometría no escapa de esta característica, de forma que, al efectuar un censo poblacional, se puede estimar un modelo definido por la Función de Regresión Poblacional (FRP), expresado de la siguiente manera:

Y_i = \beta_1 + \beta_2 X_i + u_i

Sin embargo, llevar a cabo un censo puede resultar costoso en todos los aspectos, es por esto que se recurre a muestras poblacionales, a partir de las cuales se puede estimar un modelo definido por la Función de Regresión Muestral (FRM), expresado de la siguiente manera:

\hat{Y}_i = \hat{\beta}_1 + \hat{\beta}_2 x_i + \hat{u}_i

Y si bien el objetivo principal del análisis de regresión es estimar la FRP con base en la FRM, siempre se debe tomar en cuenta que: debido a fluctuaciones muestrales, la estimación de la FRP basada en la FRM es, en el mejor de los casos, una aproximación.

Mínimos Cuadrados Ordinarios (MCO)

Entonces, los valores de \beta_1 y \beta_2 se pueden estimar a partir de una muestra, usando un modelo que cuente con el término de error \hat{u}_i más pequeño posible, sin embargo, no podemos permitir que estos errores se anulen.

El Método de los Mínimos Cuadrados Ordinarios (MCO) que consiste en considerar, de todos los modelos posibles, el modelo tal que la suma de los cuadrados de los residuos \hat{u}_i sea la más pequeña, es decir, tal que la siguiente suma sea la más pequeña:

\sum \hat{u}^2

Llevando a cabo los cálculos necesarios para que esto se cumpla, se determina que los valores que estiman a \beta_1 y \beta_2, es decir, los estimadores \hat{\beta}_1 y \hat{\beta}_2 se calculan de la siguiente forma:

El valor \beta_2 se conoce como la pendiente y su estimador es:

\hat{\beta}_2 = \dfrac{\sum x_i y_i}{\sum x_i^2}

El valor \beta_1 se conoce como el intercepto y su estimador es:

\hat{\beta}_1 = \overline{Y} - \hat{\beta}_2 \overline{X}

Conociendo estas dos expresiones, podemos hacer los cálculos correspondientes en R, pues calculando la media de las variables X y Y, usamos la siguiente sintaxis para calcular los estimadores

m.X <- mean(X)
m.Y <- mean(Y)
beta2 <- sum((X - m.X)*(Y - m.Y))/sum((X - m.X)^2)
beta1 <- m.Y - beta2*m.X

Ejemplo

Consideremos un pequeño conjunto de datos, particularmente, los datos que se encuentran en la Tabla 3.2 del libro de Econometría de Damodar N. Gujarati and Dawn Porter en su quinta edición. Este conjunto de datos proporciona los datos primarios que se necesitan para estimar el efecto cuantitativo de la escolaridad en los salarios:

ObservaciónSalarioEscolaridad
14.45676
25.777
35.97878
47.33179
57.318210
66.584411
77.818212
87.835113
911.022314
1010.673815
1110.836116
1213.61517
1313.53118
Tabla 3.2

Si queremos definir un modelo que describa el salario de una persona en función del nivel de estudio que esta persona tenga, empezamos por definir las variables salario y escolaridad, usando las siguientes instrucciones:

escolaridad <- c(6,7,8,9,10,11,12,13,14,15,16,17,18)
salario <- c(4.4567,5.77,5.9787,7.3317,7.3182,6.5844,7.8182,7.8351,11.0223,10.6738,10.8361,13.615,13.531)

Una vez definidas estas variables, podemos definir nuevas variables para almacenar la media de cada una de ellas:

m.escolaridad <- mean(escolaridad)
m.salario <- mean(salario)

Posteriormente, calculamos los estimadores:

beta2 <- sum( (escolaridad-m.escolaridad)*(salario-m.salario) )/sum( (escolaridad-m.escolaridad)^2 )
beta1 <- m.salario - beta2*m.escolaridad

Al ejecutar estas instrucciones definimos las variables y podemos ver los valores que cada una de ellas tienen, particularmente la de los estimadores que son las que nos interesan.

Modelo Lineal de Salario en función de Escolaridad | totumat.com

Mi recomendación es usar el símbolo de numeral «#» para hacer comentarios en el script y mantener orden en las instrucciones que escribimos o entender porqué las escribimos, les comparto como haría yo estas anotaciones.

Modelo Lineal de Salario en función de Escolaridad | totumat.com

Habiendo calculado los valores de los estimadores, concluimos que el modelo lineal determinado por el método de los Mínimos Cuadrados Ordinarios, que estima el comportamiento de los valores expuestos en la Tabla 3.2 es el siguiente:

Salario = -0.0144 + 0.7240 \cdot Escolaridad

Anuncios

La instrucción lm

También podemos recurrir a la instrucción lm para definir un modelo lineal, de forma que si queremos definir a la variable Y en función de la variable X, entonces usamos la siguiente sintaxis:

lm(Y ~ X)

Note que se ha usado la virguilla (~) para definir la relación entre las dos variables. Entonces, continuando con nuestro ejemplo, podemos definir el modelo lineal que describe el Salario en función de la Escolaridad usando la siguiente sintaxis:

lm(salario ~ escolaridad)

Al ejecutar esta instrucción, en la consola deberá aparecer lo siguiente:

> lm(salario ~ escolaridad)
Call:
lm(formula = salario ~ escolaridad)
Coefficients:
(Intercept)  escolaridad  
   -0.01445      0.72410  

En su pantalla debería aparecer:

Modelo Lineal de Salario en función de Escolaridad | totumat.com

MCRL | totumat.com

R: Formulario del Modelo Clásico de Regresión Lineal

A continuación se presentan las fórmulas usadas para desarrollar la teoría del Modelo Clásico de Regresión Lineal (MCRL), acompañadas por su sintaxis en R de forma respectiva.

También pudiera interesarte

Anuncios

Medias

Media de la variable X:

\dfrac{\sum X_i}{n}

m.X <- sum(X)/length(X)
# o también
m.X <- mean(X)

Media de la variable Y:

\dfrac{\sum Y_i}{n}

m.Y <- sum(Y)/length(Y)
# o también
m.Y <- mean(Y)

Estimadores

Pendiente

\hat{\beta}_2 = \dfrac{\sum x_i y_i}{\sum x_i^2}

beta2 <- sum((X - m.X)*(Y - m.Y))/sum((X - m.X)^2)

Intercepto

\hat{\beta}_1 = \overline{Y} - \hat{\beta}_2 \overline{X}

beta1 <- m.Y - beta2*m.X

Residuos

\hat{u}_i = Y_i - \hat{Y}_i

res <- Y - Y.e

\sum \hat{u}_i^2 = \sum( Y_i - \hat{\beta}_1 - \hat{\beta}_2 X_i)^2

SCR <- sum(res^2)

\sum \hat{u}_i^2 = \sum y_i^2 - \dfrac{(\sum x_i y_i)^2}{\sum x_i^2}

SCR <- sum((Y-Y.e)^2) - (sum((X-m.X)*(Y-m.Y)))^2/sum((X-m.X)^2)
Anuncios

Varianza y Error Estándar

var(\hat{\beta}_2) = \dfrac{\sigma^2}{\sum x_i^2}

ee(\hat{\beta}_2) = \dfrac{\sigma}{ \sqrt{\sum x_i^2} }

var.beta1 <- sigma2.e*sum( X^2 )/(length(X) * sum( (X-m.X)^2 ))
ee.beta1 <- sqrt(var.beta1)

var(\hat{\beta}_1) = \dfrac{ \sum X_i^2 }{n \sum x_i^2} \cdot \sigma^2

ee(\hat{\beta}_1) = \sqrt{ \dfrac{ \sum X_i^2 }{n \sum x_i^2} } \cdot \sigma

var.beta1 <- sigma2.e*sum( X^2 )/(length(X) * sum( (X-m.X)^2 ))
ee.beta1 <- sqrt(var.beta1)

\hat{\sigma}^2 = \dfrac{\sum \hat{u}_i^2}{n-2}

sigma2.e <- SCR/(lenght(X)-2)

\hat{\sigma} = \sqrt{\dfrac{\sum \hat{u}_i^2}{n-2}}

ee.e <- sqrt(sigma2.e)
Anuncios

Covarianza

cov(\hat{\beta}_1,\hat{\beta}_2) = - \overline{X} var(\hat{\beta}_2)

cov.b1b2 <- -m.X (sigma2.e/sum((X-m.X)^2))

Bondad de Ajuste

r^2 = \dfrac{SCE}{SCT} = \dfrac{\sum (\hat{y}_i - \overline{Y})^2}{\sum (Y_i - \overline{Y})^2}

r2 <- sum((Y.e - m.Y)^2)/sum((Y - m.Y)^2)

r^2 = 1 - \dfrac{SCR}{SCT} = 1 - \dfrac{ \sum \hat{u}_i^2}{\sum (Y_i - \overline{Y})^2}

r2 <- 1 - sum((Y - Y.e)^2)/sum((Y - m.Y)^2)

r^2 = \dfrac{(\sum x_i y_i)^2}{\sum x_i^2 \sum y_i^2}

r2 <- (sum((X-m.X)*(Y-m.Y))^2/(sum((X - m.X)^2)*sum((Y - m.Y)^2))

r = \dfrac{\sum x_i y_i}{\sqrt{\sum x_i^2 \sum y_i^2}}

r <- sum((X-m.X)*(Y-m.Y))/sqrt(sum((X - m.X)^2)*sum((Y - m.Y)^2))

r^2 = \dfrac{(\sum y_i \hat{y}_i)^2}{(\sum y_i^2) (\sum \hat{y}^2)}

r2 <- (sum(Y-m.Y)*(Y-Y.e))^2/(sum((Y - m.Y)^2)*sum((Y - Y.e)^2))
Anuncios

Estadísticos

Estadístico F

F = \dfrac{\hat{\beta}_2^2 \sum x_i^2}{ \hat{\sigma}^2 }

F.c <- beta2^2 * sum((X-m.X)^2)/sigma2.e
#p-value de F.c
pf(F.c,1,length(X)-2,lower.tail = F)

Estadístico t

# Dos colas
qt(alpha/2,df=length(X)-2,lower.tail = FALSE)

# Una cola
qt(alpha,df=length(X)-2,lower.tail = FALSE)

Estadístico chi cuadrado \mathbf{\chi}^2

# Cola izquierda
qchisq(alpha/2,df=length(X)-2,lower.tail = FALSE)

# Cola derecha
qchisq(1-alpha/2,df=df=length(X)-2,lower.tail = FALSE)
Anuncios

Intervalos de Confianza

\hat{\beta}_2 \pm t_{\alpha/2} ee(\hat{\beta}_2)

li.beta2 <- beta2 - qt(alpha/2, df=length(escolaridad)-2,lower.tail=F)*ee.beta2
ls.beta2 <- beta2 + qt(alpha/2, df=length(escolaridad)-2,lower.tail=F)*ee.beta2

\hat{\beta}_1 \pm t_{\alpha/2} ee(\hat{\beta}_1)

li.beta1 <- beta1 - qt(alpha/2, df=length(escolaridad)-2,lower.tail=F)*ee.beta1
ls.beta1 <- beta1 + qt(alpha/2, df=length(escolaridad)-2,lower.tail=F)*ee.beta1

P \left[ (n-2) \dfrac{\hat{\sigma}^2}{\chi^2_{\alpha/2}} \leq \sigma^2 \leq (n-2) \dfrac{\hat{\sigma}^2}{\chi^2_{1-\alpha/2}} \right] = 1-\alpha

li.sigma2 <- (length(escolaridad)-2)*sigma2.e/qchisq(alpha/2,df=length(escolaridad)-2,lower.tail=F)
ls.sigma2 <- (length(escolaridad)-2)*sigma2.e/qchisq(1-alpha/2,df=length(escolaridad)-2,lower.tail=F)
Anuncios

Predicción

\hat{Y}_o = \hat{\beta}_1 + \hat{\beta}_2 X_0

latex Y.0 <- beta1 + beta2*X.0

Varianza de la Predicción Media

var(\hat{Y}_0) = \sigma^2 \left[ \dfrac{1}{n} + \dfrac{(X_0 - \overline{X})^2}{\sum x_i^2} \right]

ee(\hat{Y}_0) = \sigma\sqrt{\dfrac{1}{n} + \dfrac{(X_0 - \overline{X})^2}{\sum x_i^2} }

varm.Y.0 <- sigma2.e*(1/length(X)+(X.0-m.X)^2/sum((X-m.X)^2))
eem.C.0 <- sqrt(varm.C.0)

Intervalo de Confianza de la Predicción Media

P \big[ \hat{\beta}_1 + \hat{\beta}_2 X_0 - t_{\alpha/2} ee(\hat{Y}_0) \leq \beta_1 + \beta_2 X_0 \leq \hat{\beta}_1 + \hat{\beta}_2 X_0 + t_{\alpha/2} ee(\hat{Y}_0) \big] = 1 - \alpha

li.C.0 <- beta1 + beta2*X.0 - qt(0.025,df=length(X)-2,lower.tail = FALSE)*eem.C.0
ls.C.0 <- beta1 + beta2*X.0 + qt(0.025,df=length(X)-2,lower.tail = FALSE)*eem.C.0

Varianza de la Predicción Individual

var(Y_0 - \hat{Y}_0) = E[Y_0 - \hat{Y}_0]^2 = \sigma^2 \left[ 1 + \dfrac{1}{n} + \dfrac{(X_0 - \overline{X})^2}{\sum x_i^2} \right]

ee(\hat{Y}_0) = \sigma\sqrt{1+\dfrac{1}{n} + \dfrac{(X_0 - \overline{X})^2}{\sum x_i^2} }

vari.Y.0 <- sigma2.e*(1+1/length(X)+(X.0-m.X)^2/sum((X-m.X)^2))
eei.C.0 <- sqrt(vari.C.0)

Intervalo de Confianza de la Predicción Media

P \big[ \hat{\beta}_1 + \hat{\beta}_2 X_0 - t_{\alpha/2} ee(\hat{Y}_0) \leq \beta_1 + \beta_2 X_0 \leq \hat{\beta}_1 + \hat{\beta}_2 X_0 + t_{\alpha/2} ee(\hat{Y}_0) \big] = 1 - \alpha

li.C.0 <- beta1 + beta2*X - qt(0.025,df=length(X)-2,lower.tail = FALSE)*eei.C.0
ls.C.0 <- beta1 + beta2*X + qt(0.025,df=length(X)-2,lower.tail = FALSE)*eei.C.0