Equipo 1


Pregunta de investigación

¿Cuál es el efecto de los cambios laborales en modadlidad virtual en la empleabilidad de las personas con discapacidad Motriz?


Contexto de la problemática

México cuenta con unos 21 millones de personas con discapacidad, de las cuales solo un millón están en capacidad de trabajar. Entre ellas, apenas 12,000 tienen un título universitario y solo 3,600 están empleadas. Menos del 0.01% de las personas con discapacidad en edad laboral ocupan puestos de gerencia o supervisión (Inegi, 2022). El 38.5% de las personas con alguna discapacidad son económicamente activas. Dentro de este grupo, las personas con discapacidad motriz representan el 30.2% (Inegi, 2020). Dentro de este grupo, se identifica que 2,939,986 personas tienen discapacidad motriz, lo que implica dificultades para caminar, subir o bajar. De estas, 1,282,534 son hombres y 1,657,452 son mujeres. Asimismo, se registra que 1,168,098 personas tienen dificultades para bañarse, vestirse o comer, otro criterio utilizado para definir la discapacidad motriz. Dentro de este subgrupo, 540,971 son hombres y 627,127 son mujeres.

Estos dos criterios, que representan dos de cada cinco criterios utilizados por el INEGI para determinar la discapacidad, son considerados como los grupos de personas con discapacidad motriz para efectos de esta investigación. Estos datos resaltan la importancia de comprender las necesidades y desafíos específicos de las personas con discapacidad motriz en el contexto laboral y social, y la urgente necesidad de promover su inclusión y equidad en la sociedad. Estas cifras reflejan una deficiencia significativa en la inclusión laboral y social de las personas con discapacidad motriz, lo que sugiere la necesidad de abordar los obstáculos y barreras que enfrentan en el mercado laboral para garantizar su plena participación, sobretodo en el sector formal.


Comportamiento histórico de variable de interés

A continuación, graficamos el comportamiento histórico de la formalidad laboral en México. Esta variable es crucial ya que el trabajo formal suele ser más accesible para personas con discapacidad. Observar su evolución nos ayuda a tomarlo como referencia para nuestro modelo.

La tasa de

Diagrama de flujo del problema

De igual forma, se construyó el diagrama de flujo del caso:

Texto alternativo

Modelo principal

library(deSolve)
library(ggplot2)

# Se establecieron las condiciones iniciales de la variable de estado
inicial.conditions<-c(PersonasConDiscapacidadMotrizEmpleadas = 3600,
                      ninos = 31226156, 
                      jovenes = 20987285, 
                      adultos = 55007819,
                      jubilados = 10209282, 
                      PoblacionConDiscapacidadBuscandoEmpleo = 100000
                      ) 

# Vector de tiempos para la simulación
times<-seq(0,25, by = 0.5)

# Función del modelo
modelo <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {

    #Variables Auxiliares o endógenas
    ninos_discapacitados <- ninos*Prop_ninos_discapacitados
    jovenes_discapacitados <- jovenes*Prop_jovenes_discapacitados
    adultos_discapacitados <- adultos*Prop_adultos_discapacitados
    jubilados_discapacitados <- jubilados*Prop_jubilados_discapacitados
    PoblacionEnEdadParaTrabajar <- jovenes_discapacitados + adultos_discapacitados
    tasa_actividad_economica <- PersonasConDiscapacidadMotrizEmpleadas * 0.0000000001
    #Se modeló de esta manera para tomar en consideración el incremento abrupto que se dió durante la pandemia. Se cosidera la pandemia como el año 5. 
    ModalidadVirtual = ifelse(t < 5, 0, TasaCoberturaModalidadVirtual)
    WhiteCollarJobs <- WhiteCollarJobs_initial * (1 + tasa_actividad_economica) #funciontabla
    VacantesLaboralesAccesibles <- WhiteCollarJobs *  (1 - tasa_discriminacion) *  (ModalidadVirtual + TasaCoberturaInfraestructuraInclusiva)
    PoblacionEconomicamenteActiva <- PoblacionEnEdadParaTrabajar * TasaDisposicionLaboral
    TasaDespidosInjustificados <- tasa_autonoma_despidos_injustificados * (1 + tasa_discriminacion) * + (1 - tasa_actividad_economica) 
    Despidos <- PersonasConDiscapacidadMotrizEmpleadas * (TasaDespidosJustificados + TasaDespidosInjustificados)
    AbandonoLaboral <- PersonasConDiscapacidadMotrizEmpleadas * (TasaRotacionNatural + TasaInsatisfaccionLaboral) 

    #Variables de flujo (son las que modifican a las variables de estado)
    pubertad <- ninos/Tiempo_Nino
    madurez <- jovenes/Tiempo_Joven
    jubilacion <- adultos/Tiempo_Adulto 
    nacimiento <- (Tasa_Fertilidad_Adulto*(adultos/2))/Tiempo_Adulto + (Tasa_Fertilidad_Joven*(jovenes/2))/Tiempo_Adulto
    muerte <- jubilados/Tiempo_Jubilado
    muertes_ninos <- ninos 
    muertes_jovenes <- jovenes 
    muertes_adultos <- adultos 
    ContratacionLaboralPersonasConDiscapacidad <- min(VacantesLaboralesAccesibles, PoblacionConDiscapacidadBuscandoEmpleo)
    DesempleoPersonasConDiscapacidad <- AbandonoLaboral + Despidos
    poblacion <- ninos + jovenes + adultos + jubilados
    Desempleados_buscan_empleo = DesempleoPersonasConDiscapacidad * Proporcion_desempleados_buscando_empleo
      
    #variable de estado
    dninos <- nacimiento-pubertad
    djovenes <- pubertad-madurez
    dadultos <- madurez-jubilacion
    djubilados <- jubilacion-muerte
    dPersonasConDiscapacidadMotrizEmpleadas <- ContratacionLaboralPersonasConDiscapacidad - DesempleoPersonasConDiscapacidad
    dPoblacionConDiscapacidadBuscandoEmpleo <- Desempleados_buscan_empleo + PoblacionEconomicamenteActiva
    
    #Devuelve los resultados de la variable de estado
    return(list(c(dPersonasConDiscapacidadMotrizEmpleadas, dninos, djovenes, dadultos, djubilados, dPoblacionConDiscapacidadBuscandoEmpleo),
                 ContratacionLaboralPersonasConDiscapacidad = ContratacionLaboralPersonasConDiscapacidad,
             PoblacionEconomicamenteActiva = PoblacionEconomicamenteActiva,
 VacantesLaboralesAccesibles = VacantesLaboralesAccesibles, 
 DesempleoPersonasConDiscapacidad = DesempleoPersonasConDiscapacidad, 
 AbandonoLaboral = AbandonoLaboral, 
 Despidos = Despidos, 
 TasaDespidosInjustificados = TasaDespidosInjustificados,
 tasa_actividad_economica = tasa_actividad_economica, 
 WhiteCollarJobs = WhiteCollarJobs, 
 nacimiento = nacimiento, 
 PoblacionEnEdadParaTrabajar = PoblacionEnEdadParaTrabajar
)) 
  })
}

# Se seleccionó el método de integración a utilizar en la simulación, en este caso 'rk4' (Runge-Kutta de 4to orden)
intg.method<-c("rk4")

# Parámetros del modelo
# Los datos utilizados están basados en el INEGI. 
# INEGI (2022) Población con discapacidad, con limitación en la actividad cotidiana 2022. Censo Poblacional 2022, INEGI.

  parameters<-c(TasaCoberturaInfraestructuraInclusiva = 20/100, 
              TasaCoberturaModalidadVirtual = 10/100,
              TasaDisposicionLaboral = 7.6/100,
              TasaDespidosJustificados = 5/100,
              TasaRotacionNatural = 5/100, 
              TasaInsatisfaccionLaboral = 5/100, 
              Tasa_Fertilidad_Adulto = 2.25, 
              Tasa_Fertilidad_Joven = 0.6333333333333, 
              Tiempo_Nino = 14, 
              Tiempo_Joven = 15, 
              Tiempo_Adulto = 31, 
              Tiempo_Jubilado = 25, 
              Prop_ninos_discapacitados = 0.01067, 
              Prop_jovenes_discapacitados = 0.01187, 
              Prop_adultos_discapacitados = 0.0293,
              Prop_jubilados_discapacitados = 0.19662, 
              tasa_autonoma_despidos_injustificados = 0.05,
              tasa_discriminacion = 0.8,
              condiciones_laborales = 0.4,
              cuotas_inclusividad = 8, 
              WhiteCollarJobs_initial = 25000000,
              Proporcion_desempleados_buscando_empleo = 0.5)
  
# Simulación utilizando la función 'ode' del paquete deSolve
out <- ode(
  y = inicial.conditions , #condiciones iniciales
  times = times, #tiempo de simulación
  func = modelo, #función del modelo
  parms = parameters ,
  method = intg.method
)

out<-as.data.frame(out)

library("ggplot2")

# Graficar los resultados de la simulación
ggplot(out, aes(x = times, y = PersonasConDiscapacidadMotrizEmpleadas)) + 
  geom_line(color = "#064F6F") +
  theme_minimal() +  # This sets a clean theme with white background and light gray grid lines
  theme(
    panel.grid.major = element_line(color = "gray80"),  # Adjust the color of major grid lines
    panel.grid.minor = element_line(color = "gray90")   # Adjust the color of minor grid lines
  )

Interpretación del resultado y gráfica:

Como se puede observar en los datos, la población con discapacidad motriz empleada experimenta un crecimiento rápido y sostenido en los primeros años, antes de estabilizarse en una tasa más moderada.

Particularmente notable es el aumento abrupto que se observa en el año 5, coincidiendo con la pandemia de COVID-19, en la cual se modeló la adaptación acelerada de las empresas al teletrabajo en respuesta a la crisis sanitaria.

A lo largo de un período de 25 años, estas dinámicas de crecimiento han llevado a un aumento significativo en la cantidad de personas con discapacidad motriz empleadas, alcanzando un total de 6 millones.

Análisis de Incertidumbres

Realizaremos varias simulaciones considerando que la tasa de infraestructura inclusiva, la adopción de modalidades virtuales y la discriminación son variables que pueden variar considerablemente. Por lo tanto, las modelarmos con diferentes valores para entender cómo influyen en la inclusión laboral de personas con discapacidad motriz.

Estas simulaciones nos permitien entender el impacto directo de estos factores en las oportunidades laborales y el ambiente de trabajo para personas con discapacidades motrices.

#=========================================================================
# Decision Analysis  
#=========================================================================
#definimos las incertidumbres  

 x1<- seq(2/100, 30/100, 1/200) #Infraestructura inclusiva. Modelamos que sea del 2% al 30%. 
 x2<- seq(2/100, 30/100, 5/200) #Modalidad virtual. Modelamos que sea del 2% al 30%.
 x3<- seq(20/100, 50/100, 5/200) #Discriminación. Modelamos que sea del 20% al 50%.

#combine all stressors 
 Xs<-expand.grid(list(TasaCoberturaInfraestructuraInclusiva=x1,TasaCoberturaModalidadVirtual=x2, tasa_discriminacion=x3))
 Xs$Run.ID <- 1:nrow(Xs)
#map all Xs into the model 
out_all <- list()
#Establecemos el loop
for (i in 1:nrow(Xs))
{
  
#Parámetros  
parameters.Xs<- parameters<-c(
              TasaCoberturaInfraestructuraInclusiva = Xs$TasaCoberturaInfraestructuraInclusiva[i], 
              TasaCoberturaModalidadVirtual = Xs$TasaCoberturaModalidadVirtual[i],
               TasaDisposicionLaboral = 7.6/100,
              TasaDespidosJustificados = 5/100,
              TasaRotacionNatural = 5/100, 
              TasaInsatisfaccionLaboral = 5/100, 
              Tasa_Fertilidad_Adulto = 2.25,
              Tasa_Fertilidad_Joven = 0.6333333333333, 
              Tiempo_Nino = 14, 
              Tiempo_Joven = 15, 
              Tiempo_Adulto = 31, 
              Tiempo_Jubilado = 25, 
              Prop_ninos_discapacitados = 0.01067, 
              Prop_jovenes_discapacitados = 0.01187, 
              Prop_adultos_discapacitados = 0.0293,
              Prop_jubilados_discapacitados = 0.19662, 
              tasa_autonoma_despidos_injustificados = 0.05,
              tasa_discriminacion = Xs$tasa_discriminacion[i],
              condiciones_laborales = 0.4,
              cuotas_inclusividad = 8, 
              WhiteCollarJobs_initial = 25000000, 
              Proporcion_desempleados_buscando_empleo = 0.15)

#Simulación del modelo
out <- ode(y = inicial.conditions,
           times = times,
           func = modelo, #recuerda actualizar el nombre de tu función para cada caso que resulevas
           parms = parameters.Xs,
           method =intg.method )

# Convertir el objeto out (que es una matriz) a un data frame para facilitar el manejo.
# Esto también permite agregar una nueva columna Run.ID en el siguiente paso.
out <- data.frame(out)

# Agregar una columna Run.ID al data frame out.
# Esto proporciona un identificador único para cada simulación.
out$Run.ID <- Xs$Run.ID[i]

# Agregar el data frame out a la lista out_all.
# out_all es una lista que recopila los resultados de todas las simulaciones.
# La función append() se utiliza para agregar el data frame out a la lista out_all.
out_all <- append(out_all, list(out))

# Imprimir el Run.ID para esta simulación particular.
# Esto podría ser útil para el seguimiento del progreso de la simulación, especialmente si hay muchas iteraciones.
# print(Xs$Run.ID[i]) POR SI QUIERES VER EN QUE SIMULACIÓN VA

} 

# Concatenar resultados
# La función do.call() ejecuta una función (en este caso "rbind") en una lista de argumentos (en este caso out_all).
# "rbind" es una función que combina data frames por filas.
# Al utilizar do.call con "rbind", se están combinando todos los data frames en la lista out_all en un solo data frame.
out_all <- do.call("rbind", out_all)

# La función dim() devuelve las dimensiones del objeto, mostrando el número de filas y columnas del data frame out_all.
# Esto puede ser útil para verificar que todos los data frames se hayan combinado correctamente.
dim(out_all)
## [1] 453492     19
# Combinar out_all y Xs
# La función merge() combina dos data frames en uno, basándose en una o más columnas en común.
# En este caso, se están combinando los data frames out_all y Xs basándose en la columna "Run.ID".
# Esto agregará las columnas de Xs al data frame out_all, alineando las filas basadas en los valores de "Run.ID".
out_all <- merge(out_all, Xs, by="Run.ID")

# Nuevamente, se utiliza dim() para mostrar las dimensiones del objeto resultante.
# Esto puede ser útil para verificar que la fusión se haya realizado correctamente y para entender cómo ha cambiado
# el tamaño del data frame como resultado de la fusión.
dim(out_all)
## [1] 453492     22
#Plot

ggplot(out_all, aes(x = time, y = PersonasConDiscapacidadMotrizEmpleadas, group = Run.ID, colour = tasa_discriminacion)) +
  geom_line() +
  scale_color_gradient(low = "blue", high = "orange") +
  theme_minimal() +
  theme(
    panel.grid.major = element_line(color = "gray80"),
    panel.grid.minor = element_line(color = "gray90")
  )

ggplot(out_all, aes(x = time, y = PersonasConDiscapacidadMotrizEmpleadas, group = Run.ID, colour = TasaCoberturaModalidadVirtual)) +
  geom_line() +
  scale_color_gradient(low = "blue", high = "orange") +
  theme_minimal() +
  theme(
    panel.grid.major = element_line(color = "gray80"),
    panel.grid.minor = element_line(color = "gray90")
  )

ggplot(out_all, aes(x = time, y = PersonasConDiscapacidadMotrizEmpleadas, group = Run.ID, colour = TasaCoberturaInfraestructuraInclusiva)) +
  geom_line() +
  scale_color_gradient(low = "blue", high = "orange") +
  theme_minimal() +
  theme(
    panel.grid.major = element_line(color = "gray80"),
    panel.grid.minor = element_line(color = "gray90")
  )

Interpretación de los resultados y gráficas:
  • Discriminación: Se observa que a mayor tasa de discriminación, menor es la empleabilidad. Por ejemplo, al reducir la tasa de discriminación del 50% al 20%, se puede observar un aumento de 1 millón de personas empleadas.
  • Tasa Cobertura modalidad virtual: La adopción de modalidades virtuales tiene un impacto significativo en la empleabilidad. Un aumento en la tasa de cobertura de modalidad virtual puede resultar en 1.5 millones de personas adicionales empleadas.
  • Tasa Cobertura Infraestructura Inclusiva. Esta variable muestra el mayor impacto. Al aumentar la tasa de cobertura de infraestructura inclusiva del 10% al 30%, la empleabilidad puede duplicarse.

En el modelo original presentado anteriormente, el cual solo era una simulación, la variable se estabiliza. Sin embargo, en estos escenarios de incertidumbre, podemos observar un crecimiento exponencial bajo ciertos valores. Esto resalta cómo las variaciones en las tasas de discriminación, cobertura de modalidad virtual y cobertura de infraestructura inclusiva pueden provocar un cambio significativo en la dinámica de la variable.

Modelo con políticas

Ahora, modelamos con el impacto que tendrían las tres políticas que estamos proponiendo:

1. Aumentar la cobertura de modalidad virtual:

Se busca incrementar la accesibilidad y las oportunidades laborales para las personas con discapacidad motriz mediante el trabajo remoto.

Estrategias:

  • Incentivos Fiscales: Ofrecer incentivos fiscales a las empresas que implementen políticas de teletrabajo inclusivas y contraten a personas con discapacidad motriz.
  • Formación y Equipamiento: Proveer programas de formación en habilidades digitales y tecnologías de la información específicas para personas con discapacidad motriz, así como subsidios para la adquisición de equipos y software adaptativo.
  • Plataformas de Empleo Virtual: Desarrollar y mantener plataformas de empleo virtual que faciliten la conexión entre empleadores inclusivos y personas con discapacidad motriz, ofreciendo soporte técnico y capacitación continua.


2. Disminuir el ISN:

Con el objetivo de incentivar a los empleadores a que contraten a personas con discapacidad motríz en sus empresas.

Estrategias:

  • Reducir el Impuesto Sobre la Nómina (ISN) de 3% a 1.5% para los salarios de las personas con discapacidad motríz, guiando así a que los empresarios maximicen sus utilidades mediante la contratación de estos. -Garantizar un espacio de “piso parejo” a las personas con discapacidad motríz, quienes muchas veces son discriminados por empleadores por sus condiciones físicas.


3. Aumentar la compensación del finiquito de 3 meses de salario a 5 meses de salario.

Con esto, se busca reducir la cantidad de despidos injustificados para personas con discapacidad motríz.

Estrategias:

  • Aumentar la cantidad de meses de sueldo en el finiquito de las personas con discapacidad motríz que son despedidas de manera injustificada, con la finalidad de evitar que sean despedidas únicamente por sus condiciones físicas. -Otorgar las garantías legales para que empleados con discapacidad motríz, a quienes se les dificulta más conseguir un empleo, puedan sobrellevar su periodo de desempleo con mayor facilidad al contar con un finiquito mayor.

Esta política asegura que las personas con discapacidad motriz cuenten con un “colchon” o garantía que les permita sobrevivir mientras encuentran otro empleo y aumente la proporción de personas en desempleo que no desistan en querer seguir buscando trabajo.


Enfoque Integral y Gradual

Estas políticas no deben ser vistas como esfuerzos aislados, sino como componentes de una estrategia integral y gradual que aborde la inclusión laboral de las personas con discapacidad motriz desde múltiples frentes:

  • Fase 1: Implementación del Teletrabajo y la Virtualidad. Este primer paso se enfoca en soluciones inmediatas y flexibles, facilitando la inserción laboral desde cualquier lugar, eliminando barreras físicas y reduciendo costos para los empleadores.

  • Fase 2: Incentivos Fiscales Directos y Reducción del ISN. Una vez establecida una base de empleo virtual, los incentivos fiscales adicionales y la reducción del ISN motivan a las empresas a contratar más personas con discapacidad motriz, promoviendo una mayor inclusión en el mercado laboral.

  • Fase 3: Protección y Garantías Legales. Finalmente, la mayor compensación en el finiquito y las garantías legales aseguran la estabilidad laboral y financiera de las personas con discapacidad motriz, creando un entorno laboral más justo y equitativo.

library(deSolve)
library(ggplot2)

# Establecer las condiciones iniciales de la variable de estado (IAN)
inicial.conditions<-c(PersonasConDiscapacidadMotrizEmpleadas = 889904,
                      ninos = 31226156, 
                      jovenes = 20987285, 
                      adultos = 55007819,
                      jubilados = 10209282, 
                      PoblacionConDiscapacidadBuscandoEmpleo = 100000
                      ) 

# Definir el vector de tiempos para la simulación
times<-seq(0,25, by = 0.5)

# Definir la función del modelo
modelo <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {
        #Aqui modelamos el efecto que tendrían nuestras políticas dentro de nuestras variables
        # Política 1: Aumentar la cobertura de modalidad virtual. Esta política tiene un efecto directo sobre nuestra variable de modalidad virtual
    efecto_politica_1 <- ifelse(t<5, 0,
                                approx(c(0, 0.20, 0.40), #Variable dependiente x = prop de empresas con modalidad virtual
                                       c(0.12, 0.30, 0.60), #Variable independiente y = efecto
                                       xout = empresas_modalidad_virtual)$y)

    # Política 2: Disminuir el ISN. Esta política tiene un efecto directo sobre nuestra variable de ContratacionLaboralPersonasConDiscapacidad
    efecto_politica_2 <- approx(c(3/100, 2.5/100, 2/100), #Variable dependiente x = Tasa ISN
                            c(0, 0.05, 0.15), #Variable independiente y = efecto
                            xout = tasa_ISN)$y

    
    # Política 3: Aumentar la compensación del finiquito. Esta política tiene un efecto directo sobre nuestra variable de TasaDespidosInjustificados
    efecto_politica_3  <- approx(c(3, 4, 5), #Variable dependiente x = meses de finiquito
                            c(0, 0.15, 0.35), #Variable independiente y = efecto
                            xout = finiquito)$y

    #Variables Auxiliares o endógenas
    ninos_discapacitados <- ninos*Prop_ninos_discapacitados
    jovenes_discapacitados <- jovenes*Prop_jovenes_discapacitados
    adultos_discapacitados <- adultos*Prop_adultos_discapacitados
    jubilados_discapacitados <- jubilados*Prop_jubilados_discapacitados
    PoblacionEnEdadParaTrabajar <- jovenes_discapacitados + adultos_discapacitados
    tasa_actividad_economica <- PersonasConDiscapacidadMotrizEmpleadas * 0.0000000001
    ModalidadVirtual = ifelse(t < 5, 0, TasaCoberturaModalidadVirtual*efecto_politica_1)
    WhiteCollarJobs <- WhiteCollarJobs_initial * (1 + tasa_actividad_economica) #funciontabla
    VacantesLaboralesAccesibles <- WhiteCollarJobs *  (1 - tasa_discriminacion) *  (ModalidadVirtual + TasaCoberturaInfraestructuraInclusiva)
    PoblacionEconomicamenteActiva <- PoblacionEnEdadParaTrabajar * TasaDisposicionLaboral
    TasaDespidosInjustificados <- tasa_autonoma_despidos_injustificados * (1 + tasa_discriminacion) * + (1 - tasa_actividad_economica) * efecto_politica_3 ## checar victor
    Despidos <- PersonasConDiscapacidadMotrizEmpleadas * (TasaDespidosJustificados + TasaDespidosInjustificados)
    AbandonoLaboral <- PersonasConDiscapacidadMotrizEmpleadas * (TasaRotacionNatural + TasaInsatisfaccionLaboral) 
    

    #Variables de flujo (son las que modifican a las variables de estado)
    pubertad <- ninos/Tiempo_Nino
    madurez <- jovenes/Tiempo_Joven
    jubilacion <- adultos/Tiempo_Adulto 
    nacimiento <- (Tasa_Fertilidad_Adulto*(adultos/2))/Tiempo_Adulto + (Tasa_Fertilidad_Joven*(jovenes/2))/Tiempo_Adulto
    muerte <- jubilados/Tiempo_Jubilado
    muertes_ninos <- ninos 
    muertes_jovenes <- jovenes 
    muertes_adultos <- adultos 
    ContratacionLaboralPersonasConDiscapacidad <- min(VacantesLaboralesAccesibles, PoblacionConDiscapacidadBuscandoEmpleo) * (1 + efecto_politica_2)
    DesempleoPersonasConDiscapacidad <- AbandonoLaboral + Despidos
    poblacion <- ninos + jovenes + adultos + jubilados
    Desempleados_buscan_empleo = DesempleoPersonasConDiscapacidad * Proporcion_desempleados_buscando_empleo
      
    #variable de estado (se establece su ecuación diferencial de ahí viene la d)
    dninos <- nacimiento-pubertad
    djovenes <- pubertad-madurez
    dadultos <- madurez-jubilacion
    djubilados <- jubilacion-muerte
    dPersonasConDiscapacidadMotrizEmpleadas <- ContratacionLaboralPersonasConDiscapacidad - DesempleoPersonasConDiscapacidad
    dPoblacionConDiscapacidadBuscandoEmpleo <- Desempleados_buscan_empleo + PoblacionEconomicamenteActiva
    
    #Devuelve los resultados de la variable de estado
    return(list(c(dPersonasConDiscapacidadMotrizEmpleadas, dninos, djovenes, dadultos, djubilados, dPoblacionConDiscapacidadBuscandoEmpleo),
                 ContratacionLaboralPersonasConDiscapacidad = ContratacionLaboralPersonasConDiscapacidad,
             PoblacionEconomicamenteActiva = PoblacionEconomicamenteActiva,
 VacantesLaboralesAccesibles = VacantesLaboralesAccesibles, 
 DesempleoPersonasConDiscapacidad = DesempleoPersonasConDiscapacidad, 
 AbandonoLaboral = AbandonoLaboral, 
 Despidos = Despidos, 
 TasaDespidosInjustificados = TasaDespidosInjustificados,
 tasa_actividad_economica = tasa_actividad_economica, 
 WhiteCollarJobs = WhiteCollarJobs, 
 nacimiento = nacimiento, 
 PoblacionEnEdadParaTrabajar = PoblacionEnEdadParaTrabajar
)) 
  })
}

# Seleccionar el método de integración a utilizar en la simulación, en este caso 'rk4' (Runge-Kutta de 4to orden)
intg.method<-c("rk4")

# Definir los parámetros del modelo
  parameters<-c(TasaCoberturaInfraestructuraInclusiva = 20/100, 
              TasaCoberturaModalidadVirtual = 10/100,
              TasaDisposicionLaboral = 7.6/100,
              TasaDespidosJustificados = 5/100,
              TasaRotacionNatural = 5/100, 
              TasaInsatisfaccionLaboral = 5/100, 
              Tasa_Fertilidad_Adulto = 2.25, 
              Tasa_Fertilidad_Joven = 0.6333333333333, 
              Tiempo_Nino = 14, 
              Tiempo_Joven = 15, 
              Tiempo_Adulto = 31, 
              Tiempo_Jubilado = 25, 
              Prop_ninos_discapacitados = 0.01067, 
              Prop_jovenes_discapacitados = 0.01187, 
              Prop_adultos_discapacitados = 0.0293,
              Prop_jubilados_discapacitados = 0.19662, 
              tasa_autonoma_despidos_injustificados = 0.05,
              tasa_discriminacion = 0.8,
              condiciones_laborales = 0.4,
              cuotas_inclusividad = 8, 
              WhiteCollarJobs_initial = 25000000,
              Proporcion_desempleados_buscando_empleo = 0.5,
              tasa_ISN = 3/100,
              finiquito = 3, 
              empresas_modalidad_virtual = 12/100)
  
# Realizar la simulación utilizando la función 'ode' del paquete deSolve
out <- ode(
  y = inicial.conditions , #condiciones iniciales
  times = times, #tiempo de simulación
  func = modelo, #función del modelo
  parms = parameters ,
  method = intg.method
)

out<-as.data.frame(out)

library("ggplot2")
#Cadena de envejecimiento
ggplot(out, aes(x=times, y= PoblacionEnEdadParaTrabajar)) + geom_line()

# Graficar los resultados de la simulación
ggplot(out, aes(x=times, y= PersonasConDiscapacidadMotrizEmpleadas)) + geom_line()

Análisis de Incertidumbres de políticas

Ahora, modelamos nuestras nuevos parámetros como incertidumbres.

#=========================================================================
# Decision Analysis  
#=========================================================================
#Definimos las incertidumbres  

 x1<- seq(12/100, 50/100, 1/100) #efecto_politica_1. Modelamos que sea de 12% a 50%
 x2<- seq(1/100, 3/100, 0.5/100) #efecto_politica_2. Modelamos que sea de 1% al 3%
 x3<- seq(3, 6, 1) #efecto_politica_3. Modelamos de 1 a 6 

#combine all stressors 
 Xs<-expand.grid(list(empresas_modalidad_virtual=x1,tasa_ISN=x2, finiquito=x3))
 Xs$Run.ID <- 1:nrow(Xs)
#map all Xs into the model 
out_all <- list()
#Establecemos el loop
for (i in 1:nrow(Xs))
{
  
#Parámetros  
parameters.Xs<- parameters<-c(
              tasa_ISN = Xs$tasa_ISN[i], 
              finiquito = Xs$finiquito[i],
              empresas_modalidad_virtual = Xs$empresas_modalidad_virtual[i],
              TasaCoberturaInfraestructuraInclusiva = 20/100, 
              TasaCoberturaModalidadVirtual = 10/100,
              TasaDisposicionLaboral = 7.6/100,
              TasaDespidosJustificados = 5/100,
              TasaRotacionNatural = 5/100, 
              TasaInsatisfaccionLaboral = 5/100, 
              Tasa_Fertilidad_Adulto = 2.25, 
              Tasa_Fertilidad_Joven = 0.6333333333333, 
              Tiempo_Nino = 14, 
              Tiempo_Joven = 15, 
              Tiempo_Adulto = 31, 
              Tiempo_Jubilado = 25, 
              Prop_ninos_discapacitados = 0.01067, 
              Prop_jovenes_discapacitados = 0.01187, 
              Prop_adultos_discapacitados = 0.0293,
              Prop_jubilados_discapacitados = 0.19662, 
              tasa_autonoma_despidos_injustificados = 0.05,
              tasa_discriminacion = 0.8,
              condiciones_laborales = 0.4,
              cuotas_inclusividad = 8, 
              WhiteCollarJobs_initial = 25000000,
              Proporcion_desempleados_buscando_empleo = 0.5)

#Simulación del modelo
out <- ode(y = inicial.conditions,
           times = times,
           func = modelo, #recuerda actualizar el nombre de tu función para cada caso que resulevas
           parms = parameters.Xs,
           method =intg.method )

# Convertir el objeto out (que es una matriz) a un data frame para facilitar el manejo.
# Esto también permite agregar una nueva columna Run.ID en el siguiente paso.
out <- data.frame(out)

# Agregar una columna Run.ID al data frame out.
# Esto proporciona un identificador único para cada simulación.
out$Run.ID <- Xs$Run.ID[i]

# Agregar el data frame out a la lista out_all.
# out_all es una lista que recopila los resultados de todas las simulaciones.
# La función append() se utiliza para agregar el data frame out a la lista out_all.
out_all <- append(out_all, list(out))

# Imprimir el Run.ID para esta simulación particular.
# Esto podría ser útil para el seguimiento del progreso de la simulación, especialmente si hay muchas iteraciones.
# print(Xs$Run.ID[i]) POR SI QUIERES VER EN QUE SIMULACIÓN VA

} 

# Concatenar resultados
# La función do.call() ejecuta una función (en este caso "rbind") en una lista de argumentos (en este caso out_all).
# "rbind" es una función que combina data frames por filas.
# Al utilizar do.call con "rbind", se están combinando todos los data frames en la lista out_all en un solo data frame.
out_all <- do.call("rbind", out_all)

# La función dim() devuelve las dimensiones del objeto, mostrando el número de filas y columnas del data frame out_all.
# Esto puede ser útil para verificar que todos los data frames se hayan combinado correctamente.
dim(out_all)
## [1] 39780    19
# Combinar out_all y Xs
# La función merge() combina dos data frames en uno, basándose en una o más columnas en común.
# En este caso, se están combinando los data frames out_all y Xs basándose en la columna "Run.ID".
# Esto agregará las columnas de Xs al data frame out_all, alineando las filas basadas en los valores de "Run.ID".
out_all <- merge(out_all, Xs, by="Run.ID")

# Nuevamente, se utiliza dim() para mostrar las dimensiones del objeto resultante.
# Esto puede ser útil para verificar que la fusión se haya realizado correctamente y para entender cómo ha cambiado
# el tamaño del data frame como resultado de la fusión.
dim(out_all)
## [1] 39780    22
#Plot

ggplot(out_all,aes(x=time,y=PersonasConDiscapacidadMotrizEmpleadas, group=Run.ID, colour=tasa_ISN))+
  geom_line()+
  scale_color_gradient(low = "blue", high = "orange")
## Warning: Removed 25140 rows containing missing values (`geom_line()`).

ggplot(out_all,aes(x=time,y=PersonasConDiscapacidadMotrizEmpleadas, group=Run.ID, colour=finiquito))+
  geom_line()+
  scale_color_gradient(low = "blue", high = "orange")
## Warning: Removed 25140 rows containing missing values (`geom_line()`).

ggplot(out_all,aes(x=time,y=PersonasConDiscapacidadMotrizEmpleadas, group=Run.ID, colour=empresas_modalidad_virtual))+
  geom_line()+
  scale_color_gradient(low = "blue", high = "orange")
## Warning: Removed 25140 rows containing missing values (`geom_line()`).

Interpretación de las gráficas:

Como se puede observar en el análisis, las políticas implementadas muestran claras correlaciones con la empleabilidad de las personas con discapacidad motriz. En primer lugar, se evidencia que a medida que aumenta la tasa de ISN (índice de seguridad nacional), la empleabilidad tiende a disminuir. Este indicador sugiere que políticas de seguridad más estrictas podrían estar relacionadas con barreras adicionales para la contratación de personas con discapacidad motriz.

Por otro lado, se ha encontrado que una menor tasa de finiquito está asociada con una mayor empleabilidad. Esto indica que prácticas laborales que facilitan la retención de empleo son beneficiosas para las personas con discapacidad motriz, permitiéndoles mantenerse en el mercado laboral de manera más sostenida.

Además, las empresas que adoptan la modalidad virtual también tienen un impacto en la empleabilidad, aunque este sea menos directo en comparación con los factores anteriores. Es claro que una mayor adopción de la modalidad virtual tiene un impacto positivo en las oportunidades laborales para personas con discapacidad motriz. Esto sugiere que la flexibilidad en el lugar de trabajo puede desempeñar un papel crucial en la inclusión laboral de este grupo demográfico.

En resumen, la combinación de políticas que promuevan un entorno laboral inclusivo, junto con medidas que faciliten la estabilidad y flexibilidad en el empleo, podría significativamente mejorar las oportunidades de empleo para las personas con discapacidad motriz. Estas conclusiones subrayan la importancia de adoptar enfoques integrales que aborden tanto los aspectos normativos como las prácticas laborales para lograr una mayor equidad en el ámbito laboral.