Ecosistema

Contexto Parte 1 y 2

Instalar paquetes y llamar librerías

#install.packages("tidyverse")
library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.2     ✔ tibble    3.3.0
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.1
## ✔ purrr     1.1.0     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

#install.packages("gplots")
library(gplots)

## 
## Attaching package: 'gplots'
## 
## The following object is masked from 'package:stats':
## 
##     lowess

#install.packages("plm")
library(plm)

## 
## Attaching package: 'plm'
## 
## The following objects are masked from 'package:dplyr':
## 
##     between, lag, lead

#install.packages("DataExplorer")
library(DataExplorer)
#install.packages("forecast")
library(forecast)

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

#install.packages("lavaan")
library(lavaan)

## This is lavaan 0.6-19
## lavaan is FREE software! Please report any bugs.

#install.packages("lavaanPlot")
library(lavaanPlot)
library(readxl)
library(lmtest)

## Loading required package: zoo
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

library(dplyr)

Tema 1. Datos de panel

Importar la base de datos

#DF para parte 1 
df_hogares <- read_excel("/Users/hugoenrique/Desktop/Universidad/8vo Semestre/Generación de Escenarios/M1/ecosistema/hogares.xlsx")

#DF para parte 2 
df_pob <- read.csv("/Users/hugoenrique/Desktop/Universidad/8vo Semestre/Generación de Escenarios/M1/ecosistema/population.csv")

#DF para parte 3
df <- read.csv("/Users/hugoenrique/Desktop/Universidad/8vo Semestre/Generación de Escenarios/M1/ecosistema/ecosistema.csv")

Importar la base de datos

summary(df_hogares)

##     HogarID            Año          Miembros       Ingreso          Gasto      
##  Min.   :  1.00   Min.   :2010   Min.   :2.00   Min.   : 5039   Min.   : 3273  
##  1st Qu.: 25.75   1st Qu.:2012   1st Qu.:3.00   1st Qu.:21198   1st Qu.:14474  
##  Median : 50.50   Median :2014   Median :4.00   Median :26854   Median :20135  
##  Mean   : 50.50   Mean   :2014   Mean   :3.52   Mean   :27066   Mean   :19879  
##  3rd Qu.: 75.25   3rd Qu.:2017   3rd Qu.:4.25   3rd Qu.:32627   3rd Qu.:25113  
##  Max.   :100.00   Max.   :2019   Max.   :5.00   Max.   :50129   Max.   :34753  
##      Ahorro          Satisfacción   
##  Min.   :-23306.6   Min.   : 1.157  
##  1st Qu.:   304.7   1st Qu.: 4.796  
##  Median :  6608.9   Median : 5.667  
##  Mean   :  7186.3   Mean   : 5.697  
##  3rd Qu.: 13555.9   3rd Qu.: 6.607  
##  Max.   : 38655.5   Max.   :10.000

str(df_hogares)

## tibble [1,000 × 7] (S3: tbl_df/tbl/data.frame)
##  $ HogarID     : num [1:1000] 1 1 1 1 1 1 1 1 1 1 ...
##  $ Año         : num [1:1000] 2010 2011 2012 2013 2014 ...
##  $ Miembros    : num [1:1000] 2 2 2 2 2 2 2 2 2 2 ...
##  $ Ingreso     : num [1:1000] 41373 49362 34228 39475 45966 ...
##  $ Gasto       : num [1:1000] 23949 25511 24558 26494 30411 ...
##  $ Ahorro      : num [1:1000] 17424 23851 9669 12981 15555 ...
##  $ Satisfacción: num [1:1000] 5.54 6.46 5.48 6.42 4.12 ...

head(df_hogares)

## # A tibble: 6 × 7
##   HogarID   Año Miembros Ingreso  Gasto Ahorro Satisfacción
##     <dbl> <dbl>    <dbl>   <dbl>  <dbl>  <dbl>        <dbl>
## 1       1  2010        2  41373. 23949. 17424.         5.54
## 2       1  2011        2  49362. 25511. 23851.         6.46
## 3       1  2012        2  34228. 24558.  9669.         5.48
## 4       1  2013        2  39475. 26494. 12981.         6.42
## 5       1  2014        2  45966. 30411. 15555.         4.12
## 6       1  2015        2  41671. 25042. 16629.         6.50

# create_report(df_hogares)
plot_missing(df_hogares)

plot_histogram(df_hogares)

plot_correlation(df_hogares)

df_hogares

## # A tibble: 1,000 × 7
##    HogarID   Año Miembros Ingreso  Gasto Ahorro Satisfacción
##      <dbl> <dbl>    <dbl>   <dbl>  <dbl>  <dbl>        <dbl>
##  1       1  2010        2  41373. 23949. 17424.         5.54
##  2       1  2011        2  49362. 25511. 23851.         6.46
##  3       1  2012        2  34228. 24558.  9669.         5.48
##  4       1  2013        2  39475. 26494. 12981.         6.42
##  5       1  2014        2  45966. 30411. 15555.         4.12
##  6       1  2015        2  41671. 25042. 16629.         6.50
##  7       1  2016        2  31272. 30256.  1016.         5.08
##  8       1  2017        2  40432. 27046. 13386.         7.17
##  9       1  2018        2  37187. 25705. 11482.         6.24
## 10       1  2019        2  41977. 25374. 16603.         6.63
## # ℹ 990 more rows

Paso 1. Generar de Datos de Panel

panel_1 <- pdata.frame(df_hogares, index = c("HogarID", "Año"))

Paso 2. Revisar Heterogeneidad

plotmeans(Ingreso ~ HogarID, main= "Heterogeneidad entre hogares para el Ingreso", data=panel_1)

Modelo 1. Regresión Agrupada (pooled)

pooled <- plm(Ingreso ~ Satisfacción, data = panel_1, model="pooling")
summary(pooled)

## Pooling Model
## 
## Call:
## plm(formula = Ingreso ~ Satisfacción, data = panel_1, model = "pooling")
## 
## Balanced Panel: n = 100, T = 10, N = 1000
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -20196.53  -5106.46   -575.98   5095.02  23468.66 
## 
## Coefficients:
##              Estimate Std. Error t-value  Pr(>|t|)    
## (Intercept)  10597.75     976.80  10.850 < 2.2e-16 ***
## Satisfacción  2890.77     166.68  17.343 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    6.8145e+10
## Residual Sum of Squares: 5.2364e+10
## R-Squared:      0.23158
## Adj. R-Squared: 0.23081
## F-statistic: 300.772 on 1 and 998 DF, p-value: < 2.22e-16

Modelo 2. Efectos Fijos (wihtin)

within <- plm(Ingreso ~ Satisfacción, data = panel_1, model="within")
summary(within)

## Oneway (individual) effect Within Model
## 
## Call:
## plm(formula = Ingreso ~ Satisfacción, data = panel_1, model = "within")
## 
## Balanced Panel: n = 100, T = 10, N = 1000
## 
## Residuals:
##       Min.    1st Qu.     Median    3rd Qu.       Max. 
## -15591.951  -3123.123    -74.284   3010.168  13134.979 
## 
## Coefficients:
##              Estimate Std. Error t-value  Pr(>|t|)    
## Satisfacción  1698.14     132.73  12.794 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    2.3013e+10
## Residual Sum of Squares: 1.9469e+10
## R-Squared:      0.15403
## Adj. R-Squared: 0.05993
## F-statistic: 163.687 on 1 and 899 DF, p-value: < 2.22e-16

Método Pooled vs Modelo de Efectos Fijos

pFtest(within,pooled)

## 
##  F test for individual effects
## 
## data:  Ingreso ~ Satisfacción
## F = 15.343, df1 = 99, df2 = 899, p-value < 2.2e-16
## alternative hypothesis: significant effects

Modelo 3. Efectos Aleatorios (random) - Método Walhus

walhus <- plm(Ingreso ~ Satisfacción, data = panel_1, model="random", random.method="walhus")
summary(walhus)

## Oneway (individual) effect Random Effect Model 
##    (Wallace-Hussain's transformation)
## 
## Call:
## plm(formula = Ingreso ~ Satisfacción, data = panel_1, model = "random", 
##     random.method = "walhus")
## 
## Balanced Panel: n = 100, T = 10, N = 1000
## 
## Effects:
##                    var  std.dev share
## idiosyncratic 23574420     4855  0.45
## individual    28789336     5366  0.55
## theta: 0.7249
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -16507.33  -3220.23   -147.96   3184.91  15215.46 
## 
## Coefficients:
##              Estimate Std. Error z-value  Pr(>|z|)    
## (Intercept)  16632.69     925.15  17.978 < 2.2e-16 ***
## Satisfacción  1831.41     131.69  13.907 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    2.6429e+10
## Residual Sum of Squares: 2.2139e+10
## R-Squared:      0.16233
## Adj. R-Squared: 0.16149
## Chisq: 193.404 on 1 DF, p-value: < 2.22e-16

Modelo 3. Efectos Aleatorios (random) - Método Amemiya

amemiya <- plm(Ingreso ~ Satisfacción, data = panel_1, model="random", random.method="amemiya")
summary(amemiya)

## Oneway (individual) effect Random Effect Model 
##    (Amemiya's transformation)
## 
## Call:
## plm(formula = Ingreso ~ Satisfacción, data = panel_1, model = "random", 
##     random.method = "amemiya")
## 
## Balanced Panel: n = 100, T = 10, N = 1000
## 
## Effects:
##                    var  std.dev share
## idiosyncratic 21631698     4651 0.393
## individual    33418160     5781 0.607
## theta: 0.7534
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -16370.54  -3188.47   -210.78   3188.52  14905.18 
## 
## Coefficients:
##              Estimate Std. Error z-value  Pr(>|z|)    
## (Intercept)  16777.35     953.98  17.587 < 2.2e-16 ***
## Satisfacción  1806.01     130.63  13.825 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    2.5757e+10
## Residual Sum of Squares: 2.1617e+10
## R-Squared:      0.16074
## Adj. R-Squared: 0.1599
## Chisq: 191.14 on 1 DF, p-value: < 2.22e-16

Modelo 3. Efectos Aleatorios (random) - Método Nerlove

nerlove <- plm(Ingreso ~ Satisfacción, data = panel_1, model="random", random.method="nerlove")
summary(nerlove)

## Oneway (individual) effect Random Effect Model 
##    (Nerlove's transformation)
## 
## Call:
## plm(formula = Ingreso ~ Satisfacción, data = panel_1, model = "random", 
##     random.method = "nerlove")
## 
## Balanced Panel: n = 100, T = 10, N = 1000
## 
## Effects:
##                    var  std.dev share
## idiosyncratic 19468528     4412 0.351
## individual    35940737     5995 0.649
## theta: 0.7733
## 
## Residuals:
##      Min.   1st Qu.    Median   3rd Qu.      Max. 
## -16275.51  -3113.76   -212.49   3188.29  14690.19 
## 
## Coefficients:
##              Estimate Std. Error z-value  Pr(>|z|)    
## (Intercept)  16869.92     981.37  17.190 < 2.2e-16 ***
## Satisfacción  1789.76     129.95  13.773 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Total Sum of Squares:    2.5332e+10
## Residual Sum of Squares: 2.1286e+10
## R-Squared:      0.15972
## Adj. R-Squared: 0.15888
## Chisq: 189.701 on 1 DF, p-value: < 2.22e-16

Modelo de Efectos Fijos vs Modelo de Efectos Aleatorios

phtest(walhus, within)

## 
##  Hausman Test
## 
## data:  Ingreso ~ Satisfacción
## chisq = 64.632, df = 1, p-value = 9.03e-16
## alternative hypothesis: one model is inconsistent

Por lo tanto, nos quedamos con el Modelo de Efectos Fijos.

Tema 2. Series de Tiempo: Mapas

Generar la Serie de Tiempo

#install.packages("remotes")
#library(remotes)
#install.packages("devtools")
#install.packages("usethis", method = "binary")
#remotes::install_github("diegovalle/mxmaps")
#install.packages("sf")
library(sf)

## Linking to GEOS 3.10.2, GDAL 3.4.2, PROJ 8.2.1; sf_use_s2() is TRUE

1

## [1] 1

library(mxmaps)
library(forecast)

Generar el mapa

df2 <- df_mxstate_2020

df_mxstate_2020$value <- df2$pop #Reemplazar aquí con tus valores

mxstate_choropleth(df_mxstate_2020)

Generar la Serie de Tiempo

df4 <- df_pob %>% filter(state == "TX")
ts <- ts(data = df4$population, start = 1900, frequency = 1)

Generar el Modelo ARIMA

arima <- auto.arima(ts)
arima

## Series: ts 
## ARIMA(0,2,2) 
## 
## Coefficients:
##           ma1      ma2
##       -0.5950  -0.1798
## s.e.   0.0913   0.0951
## 
## sigma^2 = 1.031e+10:  log likelihood = -1527.14
## AIC=3060.28   AICc=3060.5   BIC=3068.6

Generar el Pronóstico

pronostico <- forecast(arima, level = 95, h=5)
pronostico

##      Point Forecast    Lo 95    Hi 95
## 2020       29398472 29199487 29597457
## 2021       29806827 29463665 30149990
## 2022       30215183 29742956 30687410
## 2023       30623538 30024100 31222977
## 2024       31031894 30303359 31760429

plot(pronostico, main="Predicción de TX")

Tema 3. Modelos de Ecuaciones Estructurales

Contexto Parte 3

El conjunto de datos es de un ecosistema.

Las variables son:

• Calidad del Suelo:
  • SPH: pH del Suelo
  • NC: Contenido de Nutrientes
  • OM: Materia Orgánica
• Calidad del Agua:
  • CL: Niveles de Contaminantes
  • DO: Oxígeno Disuelto
  • WPH: pH del Agua
• Salud del Ecosistema
  • SD: Diversidad de Especies
  • BM: Biomasa
  • EP: Productividad del Ecosistema

Generar el modelo

modelo <- '
          # Regresiones
          EcosystemHealth ~ SoilQuality + WaterQuality
          # Variables latentes
          SoilQuality =~ SPH + NC + OM
          WaterQuality =~ CL + DO + WPH
          EcosystemHealth =~ SD + BM + EP
          # Varianzas y covarianzas
          # Intercepto
          '

Generar el Diagrama

df3 <- scale(df) #Se escala cuando los valores son muy diferentes entre sí
df4 <- cfa(modelo, df3)

## Warning: lavaan->lav_object_post_check():  
##    covariance matrix of latent variables is not positive definite ; use 
##    lavInspect(fit, "cov.lv") to investigate.

summary(df4)

## lavaan 0.6-19 ended normally after 144 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of model parameters                        21
## 
##   Number of observations                           200
## 
## Model Test User Model:
##                                                       
##   Test statistic                                17.149
##   Degrees of freedom                                24
##   P-value (Chi-square)                           0.842
## 
## Parameter Estimates:
## 
##   Standard errors                             Standard
##   Information                                 Expected
##   Information saturated (h1) model          Structured
## 
## Latent Variables:
##                      Estimate  Std.Err  z-value  P(>|z|)
##   SoilQuality =~                                        
##     SPH                 1.000                           
##     NC                  2.217    1.332    1.664    0.096
##     OM                  0.167    0.402    0.414    0.679
##   WaterQuality =~                                       
##     CL                  1.000                           
##     DO                 -0.827    0.427   -1.936    0.053
##     WPH                 0.404    0.359    1.124    0.261
##   EcosystemHealth =~                                    
##     SD                  1.000                           
##     BM                 -1.899    3.995   -0.475    0.634
##     EP                 -4.224    8.093   -0.522    0.602
## 
## Regressions:
##                     Estimate  Std.Err  z-value  P(>|z|)
##   EcosystemHealth ~                                    
##     SoilQuality        0.530    1.499    0.353    0.724
##     WaterQuality       0.583    1.507    0.387    0.699
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   SoilQuality ~~                                      
##     WaterQuality     -0.079    0.051   -1.558    0.119
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .SPH               0.926    0.103    9.016    0.000
##    .NC                0.657    0.226    2.906    0.004
##    .OM                0.993    0.100    9.978    0.000
##    .CL                0.953    0.118    8.086    0.000
##    .DO                0.966    0.108    8.977    0.000
##    .WPH               0.988    0.099    9.938    0.000
##    .SD                0.992    0.100    9.954    0.000
##    .BM                0.985    0.102    9.618    0.000
##    .EP                0.948    0.166    5.720    0.000
##     SoilQuality       0.069    0.057    1.197    0.231
##     WaterQuality      0.042    0.075    0.558    0.577
##    .EcosystemHelth    0.018    0.073    0.246    0.806

lavaanPlot(df4, coef=TRUE, cov=TRUE)