Aplicación de técnicas estadísticas, matemáticas y de Inteligencia Artificial para la modelización ecológica

Aplicación de técnicas estadísticas, matemáticas y de Inteligencia Artificial para la modelización ecológica TRABAJO FIN DE GRADO Doble Grado de Matemáticas y Estadística Magdalena Jáñez Vaz Sevilla, Junio de 2023 tutorizado por Luis Valencia Cabrera Índice general Resumen........................................ iii Abstract........................................ iv ÍndicedeFiguras................................... vi ÍndicedeTablas....................................viii I Introducción 1 1. Introducción 3 1.1. Motivación.................................... 3 1.2. Objetivos .................................... 4 1.3. Estructura del documento . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 2. Preliminares 7 2.1. Modelizaciónecológica ............................. 7 2.2. Procambarusclarkii............................... 8 2.3. Técnicas de análisis descriptivo . . . . . . . . . . . . . . . . . . . . . . . . 10 2.4. Métodosdeinferencia.............................. 12 2.5. Software..................................... 18 II Estudio estadístico de los datos 23 3. Importación, tratamiento y descripción de los datos 25 3.1. Archivo OCC: ocurrencias individuales . . . . . . . . . . . . . . . . . . . . 25 3.2. Archivo MOF: información biométrica . . . . . . . . . . . . . . . . . . . . 27 3.3. ArchivoEV:eventos .............................. 28 3.4. ArchivoDatabases ............................... 30 3 4ÍNDICE GENERAL 4. Análisis descriptivo 35 4.1. Faseexploratoria ................................ 35 4.2. Diferencias entre hábitats . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 4.3. Variablesambientales.............................. 38 4.4. Tamaño según sexo y madurez . . . . . . . . . . . . . . . . . . . . . . . . . 40 4.5. Capturasentrampas.............................. 41 4.6. Proporción de individuos inmaduros . . . . . . . . . . . . . . . . . . . . . . 48 III Modelos predictivos 51 5. Introducción a los modelos predictivos 53 5.1. Ideasiniciales.................................. 53 5.2. Análisis distribucional . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 5.3. Modelos de clasificación . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 5.4. Modelosderegresión.............................. 59 6. Modelos de clasificación 63 6.1. PrediccióndelSexo............................... 63 6.2. Predicción de la Madurez . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 7. Modelos de regresión 75 7.1. Tamaño de los cangrejos . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 7.2. Número de capturas en trampas . . . . . . . . . . . . . . . . . . . . . . . . 86 IV Conclusiones 91 8. Conclusiones 93 8.1. Aportaciones iniciales . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93 8.2. Conclusiones................................... 94 8.3. Caminosdescartados.............................. 96 8.4. Líneas futuras de trabajo . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 A. Apéndice: Código omitido 99 A.1.Importación................................... 99 A.2.Análisisdescriptivo...............................102 A.3.Modelosclasificación ..............................104 A.4.Modelosregresión................................107 B. Apéndice: Código descartado 117 B.1. Regresión para capturas en trampas . . . . . . . . . . . . . . . . . . . . . . 117 B.2. Proporción de machos y hembras . . . . . . . . . . . . . . . . . . . . . . . 126 B.3. Tuning bosques aleatorios . . . . . . . . . . . . . . . . . . . . . . . . . . . 137 Bibliografía 141 i Resumen La modelización ecológica trata de comprender y traducir las complejas interacciones presentes en un ecosistema y expresarlas como un modelo matemático de forma más esquematizada de cara a ayudarnos a comprender mejor la realidad que estamos estudiando, explicar mejor su comportamiento ante determinados escenarios de interés y tratar de predecir la respuesta del sistema ante nuevos casos. El presente proyecto trabaja con datos de campo recogidos en el parque natural de Doñana, Huelva, relativos al Cangrejo Rojo Americano, Procambarus clarkii. En él se intentará procesar la información recolectada para presentarla de forma más clara y sencilla, descifrar el comportamiento de las variables implicadas, sus interacciones, vislumbrar patrones de respuesta temporal, así como utilizarla para modelizar distintos aspectos de esta especie y del ecosistema. Al tratarse de una especie invasora, es de interés analizar formas de controlar su población, y poder comprender mejor el desarrollo y comportamiento que tiene en las localizaciones estudiadas. Entre los modelos con los que se trabajará hay algunos dedicados a predecir las capturas en trampas de cangrejos según las condiciones climáticas, lo que nos sirve como una posible estimación del tamaño de la población, así como para poder observar cuál es la tendencia local, global o ante distintos posibles escenarios. Otros modelos se dedican a estudiar características de individuos en concreto como su tamaño, sexo y etapa de madurez. Por último, simplemente matizar que este estudio trata de atacar un problema real, con todo lo que ello conlleva. Al tratarse de una situación auténtica, el interés es mayor, ya que puede ayudar a comprender mejor un problema existente que se está enfrentando. Por otro lado, implica un aumento en la dificultad debido a las limitaciones de los datos frente a los conjuntos de datos que suelen darse en problemas ficticios ya preparados para realizar el estudio. Palabras clave: cangrejo rojo americano, Procambarus clarkii, Doñana, modelización ecológica. iii Abstract Ecological modeling is all about understanding and translating the complex interactions present in an ecosystem and express them as a mathematical model in a more schematic way in order to help us better understand the reality we are studying, better explain its behavior in certain scenarios of interest and try to predict the response of the system to new cases. The present project works with field data collected in the natural park of Doñana, Huelva, related to the Red Swamp Crayfish, Procambarus clarkii. In this project we will try to process the information collected to present it in a clearer and simpler way, to decipher the behavior of the variables involved, their interactions, to glimpse patterns of temporal response, and to use it to model different aspects of this species and the ecosystem. Since it is an invasive species, it is of interest to analyze ways to control its population, and to better understand its development and behavior in the locations studied. Among the models we will work with, there are some dedicated to predict the catches in crab traps according to climatic conditions, which serves as a possible estimate of the population size, as well as to observe what the local and global trend is, as well as in the face of different possible scenarios. Other models are dedicated to study characteristics of specific individuals such as their size, sex and maturity stage. Finally, it is important to note that this study attempts to address a real problem, with all that this entails. Since it is a real situation, the interest is greater, as it can help to better understand an existing problem that is being faced. On the other hand, it implies an increase in difficulty due to the limitations of the data as opposed to the data sets that usually occur in fictitious problems already prepared for the study. Key words: Red Swamp Crayfish, Procambarus clarkii, Doñana, ecological modeling. iv Índice de figuras 2.1. Tabla comparativa características de los cangrejos de interés. Fuente: [23] . 9 2.2. Foto comparativa aspecto de los cangrejos de interés. Fuente: [8] . . . . . . 10 4.1. Captura documento shiny . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 4.2. Localizaciones geográficas . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 4.3. Longitud de los cangrejos según el hábitat . . . . . . . . . . . . . . . . . . 37 4.4. Variables ambientales más correlacionadas . . . . . . . . . . . . . . . . . . 38 4.5. Evolución de la temperatura . . . . . . . . . . . . . . . . . . . . . . . . . . 39 4.6. Conductividad según hábitat . . . . . . . . . . . . . . . . . . . . . . . . . . 39 4.7. Distribución de las medidas de tamaño según sexo del cangrejo . . . . . . . 40 4.8. Distribución de las medidas de tamaño según madurez del cangrejo . . . . 41 4.9. Distribución capturas en trampas . . . . . . . . . . . . . . . . . . . . . . . 42 4.10. Distribución capturas en trampas según hábitat . . . . . . . . . . . . . . . 43 4.11.Capturasponderadas.............................. 44 4.12. Capturas en localizaciones más estudiadas . . . . . . . . . . . . . . . . . . 44 4.13. Capturas ponderadas agrupadas por estación . . . . . . . . . . . . . . . . . 45 4.14. Capturas por estación en localizaciones más estudiadas . . . . . . . . . . . 46 4.15. Capturas según estación y hábitat . . . . . . . . . . . . . . . . . . . . . . . 46 4.16. Capturas según estación y hábitat . . . . . . . . . . . . . . . . . . . . . . . 48 4.17. Proporción de inmaduros capturados por mes . . . . . . . . . . . . . . . . 49 4.18. Individuos inmaduros y totales capturados por mes . . . . . . . . . . . . . 49 4.19. Individuos inmaduros por estación . . . . . . . . . . . . . . . . . . . . . . . 50 7.1. Correlaciones variables tamaño . . . . . . . . . . . . . . . . . . . . . . . . 77 7.2. Tamaños cefalotórax según fecha . . . . . . . . . . . . . . . . . . . . . . . 77 7.3. Largo cuerpo y peso según fecha . . . . . . . . . . . . . . . . . . . . . . . . 78 7.4. Correlaciones variables tamaño según periodos temporales . . . . . . . . . 78 v 1.2. OBJETIVOS 1.2. Objetivos El principal objetivo en el estudio de una especie invasora es aprender cómo controlar su población y cómo reducir el impacto que tiene sobre el ecosistema invadido. En el caso del Cangrejo rojo americano, el tiempo que lleva asentado en Doñana hace que su población sea bastante estable, hasta el punto de que su desaparición podría implicar nuevos cambios en el ecosistema cuyas consecuencias también habría que estudiar, y analizar no solo los beneficios de su supresión sino otros nuevos problemas que llevaría emparejados. Aun así, ser capaces de estimar la cantidad de cangrejos puede ser útil de cara a paliar los efectos negativos que pueda tener, o potenciar los positivos, por lo que este es uno de los objetivos del estudio. Además, el diseño de modelos nos ayudará a entender mejor los entresijos de las relaciones entre los fenómenos involucrados, los individuos que interactúan en el ecosistema, los elementos propios del mismo, los procesos derivados de la biología de la especie y los posibles factores bióticos y abióticos implicados. Otro aspecto que sería de útilidad sería poder clasificar en función del tamaño de los cangrejos su sexo, de modo que no fuera necesario examinarlos uno a uno si se obtuviera un modelo lo suficientemente confiable. Lo mismo con la clasificación entre crías o individuos maduros, esta es otra de las variables que registraban los ecólogos, y puede ayudar a estimar el envejecimiento de la población. Para poder predecir la cantidad de huevos que se pondrán en cada puesta sería necesario poder estimar el número de hembras maduras que hubiera simultáneamente, para ello podrían utilizarse los dos modelos anteriormente planteados. De cara a la cadena alimenticia de Doñana, ya se ha comentado que el cangrejo está establecido en este ecosistema, por tanto ha pasado a formar parte de la cadena y a ser utilizado como alimento por otras especies que también habitan en el parque natural. Por ende, es de interés predecir el tamaño que alcanzarán los cangrejos, ya que eso ayuda a predecir la cantidad de alimento que ellos pueden suponer para el resto de especies. 1.3. Estructura del documento Este documento se compone de cuatro bloques distintos, cada uno compuesto a su vez de capítulos. Antes de empezar el primer bloque se encuentran el resumen, su traducción al inglés y los índices de tablas y figuras. El bloque I, Introducción, se compone de dos capítulos. Este capítulo de introducción, el capítulo 1, forma parte de él junto con el capítulo 2. Dicho capítulo trata sobre los preliminares a este estudio, en él se habla sobre modelización ecológica en general y sobre el cangrejo rojo americano, comentando tanto de su llegada a España como de aspectos propios de la especie. Se explican también las herramientas utilizadas: las técnicas de análisis descriptivo, los métodos de inferencia utilizados y el software con el que se llevaron a cabo ambas cosas. El bloque II, Estudio estadístico de los datos, contiene los capítulos 3 y 4. El capítulo 3 se dedica a la importación de los datos y su correspondiente tratamiento, se hace una breve introducción a cada conjunto tratado, así como una descripción de cada una de las variables que los componen. Una vez se han tratado los datos de forma adecuada y se han reducido a aquella información que se utilizará en el estudio, se muestra la forma en 4CAPÍTULO 1. INTRODUCCIÓN 1.3. ESTRUCTURA DEL DOCUMENTO que ha quedado dicho conjunto mediante un resumen de las variables. El capítulo 4 trata del análisis descriptivo, se hace hincapié en determinados aspectos sobre las variables, la distribución que siguen y las interacciones entre ellas. La forma de mostrarlo es principalmente mediante gráficas propiamente explicadas, aunque también se incluyen tablas numéricas para dar mayor firmeza a las conclusiones extraídas. El bloque III, Modelos predictivos, también se compone de tres capítulos. El capítulo 5 es una introducción, en la que se explican las decisiones tomadas durante el estudio realizado, para así en los dos siguientes capítulos centrarse en los modelos que finalmente fueron seleccionados. Además de explicar por qué se descartaron algunos enfoques de modelización y hacer un análisis distribucional de una de las variables, en este capítulo se describen y analizan bloques de código que se utilizarán en los capítulos siguientes. El capítulo 6 trata sobre los modelos de clasificación; para cada una de las variables categóricas se realizan varios modelos distintos, se muestran y se comparan los resultados obtenidos. El capítulo 7, por su parte, trata sobre los modelos de regresión. En este caso no es tan sencillo como aplicar el modelo a los datos, sino que requiere de una toma de decisiones que se detalla en este capítulo, así como los resultados y conclusiones obtenidas. El bloque IV, Conclusiones, se compone del capítulo 8 dedicado a las conclusiones, último del núcleo del trabajo, así como de dos apéndices. En dicho último capítulo, se resumen tanto las contribuciones realizadas como las lecciones aprendidas y principales conclusiones extraídas del estudio y desarrollo del trabajo. Adicionalmente, se analizan posibles vías de trabajo futuro. Por su parte en el primer apéndice se recopila el código utilizado a lo largo del estudio pero que no fue mostrado en el resto de la memoria, de modo que se ha destacado en el núcleo de la misma los aspectos esenciales pero se proporciona aquí en forma de apéndice por si se quiere consultar cualquier aspecto a nivel técnico. En el segundo apéndice se recopila el código que fue descartado por no dar resultados suficientemente buenos pero que merecía la pena ser mencionado. De nuevo, se proporciona como apéndice para permitir su consulta en caso de que resulte de interés de los caminos descartados o de recopilación de pruebas realizadas, pero queda fuera del foco principal del trabajo por no haber arrojado ninguna conclusión fuerte generalizable sobre caminos futuros a descartar, ya que un resultado en ese sentido también sería de interés para evitar explorar esa vía por parte de futuros trabajos por parte de investigadores y/o estudiantes. CAPÍTULO 1. INTRODUCCIÓN 5 Capítulo 2 Preliminares El capítulo a continuación está dedicado a establecer las bases teóricas sobre las que se asentará el trabajo realizado. Entre otras cosas, se habla de una posible ruta a seguir para la modelización ecológica, se introduce la especie sobre la que se estará tratando, se describen las técnicas utilizadas y el software con el que se implementaron. 2.1. Modelización ecológica La modelización ecológica es, como toda modelización, una simplificación de un sistema mucho más complejo. Los ecosistemas involucran una gran cantidad de procesos que tienen lugar, individuos que interactúan entre sí y con el entorno, y todo ello supone una gran complejidad y diversidad estructural. Además, su estudio suele abarcar grandes periodos de tiempo y espacios tan amplios que hacen difícil trabajar con el sistema sin tener que reducirlo a términos más simples. Modelizar las interacciones entre las especies que conviven en un mismo ecosistema requiere de un conocimiento sobre la conducta y comportamiento de cada una, así como sobre el clima y las condiciones del ecosistema. Blanco [4] explica los pasos a seguir en la modelización ecológica, comenzando por crear un mapa conceptual en el que se describan los elementos que componen el sistema, se identifiquen las componentes y los procesos que las unen. Una vez se tiene un esquema de los elementos que se estudiarán, el siguiente paso es construir un diagrama de flujo en el que se refleje la influencia que tienen entre sí, así como las variables reguladoras que influencian estos flujos. Es importante que este paso evitar un exceso de complejidad. El tercer paso es crear un diagrama de bucles causales, pudiendo representar así la retroalimentación entre variables del modelo (que puede ser negativa o positiva). Esto es una herramienta de comunicación más que una herramienta analítica. A continuación se estiman los parámetros y condiciones iniciales, esto es, la calibración del sistema, que debe realizarse con datos de campo. En este paso el modelo pasa de ser conceptual a ser más numérico, y por tanto requiere de una gran atención al detalle a la hora de elegir los valores o el rango que se asignará a cada variable. En el siguiente paso el modelo se ejecuta varias veces, añadiendo variaciones a los valores antes escogidos para ver la influencia que un pequeño cambio puede tener. Esto es el análisis de sensibilidad del modelo. Por último se lleva a cabo la evaluación y validación del modelo, en el que se comprueba si las predicciones del modelo son lo suficientemente cercanas a datos reales independientes del modelo. 7 2.2. PROCAMBARUS CLARKII La magnitud de esta tarea sobrepasa al alcance de este trabajo, y no es el foco principal del mismo. Además, la información disponible no es suficiente para conseguir un modelo como el explicado anteriormente. No obstante, para llevar a cabo varios de los principales pasos de interés en ese proceso descrito o en un flujo similar, es conveniente profundizar en el conocimiento de todas las “entidades” involucradas, la información disponible, y poder incrementar el conocimiento que se tiene a partir de los datos en bruto recopilados por trabajadores de campo. Esta es la idea que se planteó en este trabajo, el estudio estadístico de la información disponible, con idea de extraer conocimiento a partir de ella, que sea útil para el diseño y la validación de unos primeros modelos explicativos y predictivos basados en machine learning y para caminar hacia esos modelos ecológicos de un mayor nivel de complejidad pero que se verían altamente beneficiados por las conclusiones de nuestros estudios. De modo que se llevará a cabo una versión reducida en la que se estudiará la población del cangrejo rojo americano asentada en el parque natural de Doñana. 2.2. Procambarus clarkii El cangrejo rojo americano (Procambarus clarkii, según fue nombrado por Charles Girard en 1852) es una especie crustáceo decápodo originaria del sureste de Estados Unidos y el noreste de México. Tiene una gran capacidad de adaptación a distintos ecosistemas, tolerando incluso aguas algo salinas lo cual es inusual siendo un cangrejo de río, lo que le ha permitido extenderse por otros continentes en los que se le considera especie invasora. Junto con el cangrejo señal (Pacifastacus leniusculus), es una de las especies invasoras más representativas, habiendo llegado a prácticamente todos los continentes. España es uno de los mayores productores mundiales, siendo China el principal productor de esta especie a día de hoy [1]. Según [3] con datos hasta 2019 se importan más de 3000 toneladas al año, lo cual genera unos beneficios anuales de 3.300.000 €al año. Por otro lado, según [2], con datos hasta 2008, eran alrededor de 2 millones al año, y proporcionaba empleo para 400 personas en plantas de procesado de marisco y varios cientos de pescadores. Su llegada a España se produjo de forma premeditada con fines comerciales. Las diversas fuentes coinciden en que ocurrió a principios de la década de los 70, y aunque no llegan a ponerse de acuerdo con el año exacto se estima que fue entre 1969 y 1973. Se propuso a raíz de los buenos resultados que la cría de esta especie había tenido en Estados Unidos, generando una industria millonaria. Los expertos que asesoraron sobre el tema contaban con poco conocimiento sobre el funcionamiento de las enfermedades de los cangrejos en esa época y se consideró que el impacto ecológico que tendría la introducción de esta especie sería mínimo. Los primeros ejemplares de esta especie que llegaron a la Península fueron importados legalmente desde Monroe (Louisiana, EE.UU.) hasta cerca de Puebla del Río, en las marismas del Guadalquivir. [23] Los cangrejos no eran criados en granjas de astacicultura sino que se capturaban en su hábitat natural. Esto, junto con la sencillez para transportarlos, provocó que su expansión fuera rápida y poco controlada, ya que los pescadores o particulares que querían obtener beneficios por la venta de esta especie podían llevar a cabo introducciones en nuevas zonas. Por estos motivos el cangrejo rojo no tardó en ocupar la mayor parte del bajo Guadalquivir, de modo que hacia 1976 se comercializaban desde esta zona unos 9.000 kg de cangrejos, y en 1979 la cantidad era diez veces mayor [23]. A día de hoy los únicos lugares de la Península en las que parece escasear son Galicia y algunas otras zonas del 8CAPÍTULO 2. PRELIMINARES 2.2. PROCAMBARUS CLARKII Figura 2.1: Tabla comparativa características de los cangrejos de interés. Fuente: [23] norte. A pesar de lo que se creía en un principio, la introducción del cangrejo rojo americano en los ecosistemas acuáticos de la Península Ibérica ha provocado enormes cambios en sus características y funcionamiento. Algunos de estos efectos han sido la reducción o eliminación de vegetación acuática sumergida y transformar las aguas claras de los ecosistemas en que habitan en sistemas de aguas turbias, debido a la movilización de sedimentos del fondo que suelen llevar a cabo. Aparte de eso, el cangrejo rojo americano, al igual que el cangrejo señal, es portador de la afanomicosis o peste del cangrejo y, aunque se ve afectado por ella con baja frecuencia, los cangrejos de río europeos sí que son más susceptibles a esta enfermedad. La alimentación del cangrejo es omnívora y oportunista, pudiendo alimentarse de plantas así como de otros macroinvertebrados acuáticos, y larvas o puestas de peces y anfibios, lo que suele ser propio de individuos juveniles. Al mismo tiempo, el cangrejo ha supues- to una fuente de alimento para muchas especies de animales como pueden ser anguilas, cigüeñas, gaviotas, milanos, zorros o nutrias. Esta última, la Lutra lutra, es una especie vulnerable según el Catálogo Español de Especies Amenazadas, de la cual el cangrejo puede suponer un 76% de su alimentación si la población es suficientemente abundante. En general, las especies depredadoras del cangrejo han aumentado su riqueza en las marismas del Guadalquivir de manera considerable respecto a las especies herbívoras o que tienen un menor porcentaje de cangrejo en su dieta [25]. En cuanto a las condiciones necesarias para vivir, habitualmente los cangrejos de río son exigentes en cuanto a las características de las aguas en las que se establecen, pero en comparación el cangrejo rojo, que nos ocupa en este trabajo, tiene una gran tolerancia ecológica y flexibilidad a la hora de adaptarse a las condiciones del ecosistema. En general prefiere aguas templadas, y acostumbra a excavar refugios en los que enterrarse en periodos desfavorables como pueden ser sequías o épocas frías. Gracias a esto puede resistir en seco hasta 4 meses; también es resistente a la contaminación orgánica, las bajas concentraciones de oxígeno disuelto y a un amplio rango de salinidades. Por contra, destacar que parece ser delicado respecto a la altitud y otros factores que esta conlleva, siendo escaso por encima de 600 metros. El aspecto reproductivo es una de las principales diferencias entre el cangrejo rojo CAPÍTULO 2. PRELIMINARES 9 2.3. TÉCNICAS DE ANÁLISIS DESCRIPTIVO Figura 2.2: Foto comparativa aspecto de los cangrejos de interés. Fuente: [8] americano y otros como el cangrejo señal o el cangrejo de río europeo, Austropotamobius pallipes. Este último es el cangrejo nativo que está siendo desplazado por el cangrejo rojo americano, y actualmente es una especie amenazada. Para comparar las principales diferencias entre estos tres cangrejos, se dispone de la figura 2.1. Como puede verse en ella, el cangrejo rojo americano, P. clarkii, es mucho menos longevo que la especie nativa, A. pallipes, o que otro de los principales cangrejos de río invasores, P. leniusculus, pero a cambio tiene una maduración prácticamente inmediata: mientras que el resto tardan más de dos años en alcanzar la madurez, el cangrejo rojo americano no tarda ni medio año. El número de huevos es bastante superior al del cangrejo local, estando igualado con el cangrejo señal, pero de nuevo destaca por el poco tiempo que necesitan incubarse. Aparte de su rápido crecimiento y maduración sexual, otro factor a tener en cuenta es su capacidad de adaptar sus ciclos a condiciones climáticas y ambientales muy diversas. Si las condiciones son favorables, prácticamente en cualquier época del año pueden existir hembras con huevos, y estas pueden reproducirse más de una vez en una temporada. Combinando esto con su alta resistencia a enfermedades, no es difícil entender el motivo de que tenga tango éxito a la hora de extenderse por nuevos ecosistemas. A pesar de todos los cambios e inconvenientes que acarreó la introducción del cangrejo rojo americano en la Península Ibérica, a día de hoy lleva ya casi cinco décadas en la zona, y se ha establecido de forma bastante estable. Ha encontrado un lugar en la cadena alimenticia de los ecosistemas en los que habita, hasta el punto de que su desaparición podría acarrear problemas para las especies que se alimentan de él. Es por eso que algunos estudios han desistido de intentar disminuir su población para enfocarse en intentar mitigar los efectos que pueda tener su convivencia con aquellas especies con las que cohabitan e intentar sacar provecho de esta. 2.3. Técnicas de análisis descriptivo Una vez se tienen los datos en el formato adecuado, es importante saber mirar qué estructura tienen, estudiar el rango de valores que toma cada variables, así como la interacción que puedan tener entre sí. Aunque el objetivo final sea obtener estimaciones y modelos que ajusten las variables de interés, un estudio previo permite que eso no tenga 10 CAPÍTULO 2. PRELIMINARES 2.3. TÉCNICAS DE ANÁLISIS DESCRIPTIVO que hacerse a ciegas. Es necesario por tanto saber qué mirar, y más aún qué es mejor ignorar, para tener una visión más clara de los datos con los que se trabaja. En el capítulo de Análisis descriptivo se representan variables en distintos formatos. En esta sección se justifica la elección de cada formato en función de lo que quería resaltarse con las gráficas o tablas añadidas. Los gráficos de dispersión, o scatter plots, son útiles para enfrentar dos variables numéricas, cada una colocada en uno de los ejes. Para estudiar la correlación entre dos variables es la mejor opción, pues cada punto tiene una coordenada de cada variable, como puede verse en la Figura 4.4. Si las variables son independientes la distribución es uniforme, mientras que en caso de tener una alta correlación se agruparían los puntos a lo largo de una línea (la pendiente dependería de si la correlación es positiva o negativa). Si a esto se suma un ajuste como el que proporciona geom_smooth de la librería ggplot2, se ve de forma todavía más clara dicha recta en caso de que las variables estén altamente correladas. A los gráficos de dispresión se les pueden añadir estéticas, como es el color o la forma de los puntos en función de otra variable distinta a las que se usaron para los ejes. Generalmente se usan variables categóricas, ya que así cada color o forma corresponde a una de las categorías, y es posible ver si las que pertenecen a un mismo grupo se encuentran juntas o si están dispersas. El color proporciona un matiz más visual, es preferible sobre todo cuando se dispone de un conjunto de tamaño considerable. Sin embargo ante la posibilidad de que un documento sea consultado impreso en blanco y negro, conviene aplicar también la estética de la forma de los puntos. El tamaño de los puntos o su transparencia son también aspectos que se pueden modificar, aunque esto sirve más para cambiar la estética de todos los puntos a la vez más que para hacerlo en función de otra variable. En la Figura 4.14 se recurre a esto para poder diferenciar las múltiples líneas de la gráfica. Un caso concreto en el que los gráficos de dispersión son altamente beneficiosos es para comparar los valores reales de una variable y las predicciones que se obtienen mediante un modelo predictivo. En este caso es imprescindible que para cada individuo se pueda visualizar al mismo tiempo el valor real y el estimado, así como lo cerca que están ambos valores, por lo que utilizar cada valor para una coordenada es apropiado. El ajuste perfecto ocurriría si todos los valores representados estuvieran sobre la diagonal (del primer cuadrante en caso de trabajar con valores positivos), por lo que es conveniente al mismo tiempo que se representan los pares de valores representar también la diagonal. Esto suele hacerse con una línea discontinua, por cuestiones estéticas. Puede verse un ejemplo de esto en la Figura 7.10. Por último, si se dispone de las coordenadas de lugares de interés (como ocurre en este estudio con las localizaciones de capturas de cangrejos), al representarlas de esta forma enfrentando Latitud y Longitud puede verse la silueta del mapa sobre el que se está trabajando. Puede verse en la Figura 4.2. Cuando solo quiere estudiarse una única variable numérica, si solo toma valores enteros los diagramas de barras ofrecen una buena estimación de su función densidad. En el eje Y se realiza el conteo de la frecuencia con que la variable aparece en cada intervalo representado. Si la variable no solo toma valores enteros puede usarse también esta representación, denominada en el caso de variables continuas histograma, teniendo en cuenta que es importante la elección del número de divisiones realizado. Un ejemplo de esto se ve en la Figura 4.9. Representar la función densidad es también una buena forma de visualizar la CAPÍTULO 2. PRELIMINARES 11 2.4. MÉTODOS DE INFERENCIA distribución que siguen los valores de la variable. Para comparar una variable numérica y una variable categórica se pueden usar diagramas de caja y bigote o diagramas de violín, que son bastante similares entre sí, o crear el diagrama de barras como se mencionó anteriormente, pero generando uno para cada valor de la variable categórica. En el primer tipo de gráfico se representan los cuartiles, así como los valores atípicos de la variable numérica para cada valor de la variable categórica. En el segundo tipo de gráfico la representación es más similar a la función de densidad, teniendo una forma simétrica y más redondeada que la caja. Se pueden añadir los cuantiles deseados al diagrama de violín, por lo que la información que puede proporcionar es más completa, aunque para casos simples es preferible el de caja y bigotes, que muestra solo la información básica. Al representarse lado a lado las funciones de densidad para cada categoría es sencillo compararlas a ojo, lo que permite sopesar si es conveniente aplicar un test estadístico que compare ambas distribuciones y confirmar si son distintas, o si al contrario el grupo al que se pertenece no influye en esta variable lo suficiente y por tanto es mejor trabajar con el conjunto al completo. Ejemplos de esto son la Figura 4.3, la segunda gráfica de la Figura 4.5 y la Figura 4.15. Para enfrentar dos variables categóricas las gráficas no son de tanta utilidad, pero las tablas de contingencia permiten mostrar rápidamente el conteo de cada posible cruce entre las categorías de cada variable. También sirve por si las variables tienen valores perdidos, pues estos se pueden contar como una categoría más y se comprueba de qué forma se distribuyen los valores perdidos de cada una respecto a las categorías de la otra. Esto se utiliza en la Tabla 4.1. Las tablas con datos de resumen presentando distintos estadísticos también pueden servir para mostrar valores como la media, otros cuantiles o la desviación típica de una variable en función de varias categorías, al igual que se hacía con los diagramas de caja y bigotes o de violín, pero pudiendo ver el valor numérico. En ocasiones es suficiente solo ver la gráfica, pero en otras la información numérica acompañada del gráfico puede ser aún más revelador. La Tabla 4.2 es un ejemplo de esto. 2.4. Métodos de inferencia Tras haber realizado el análisis descriptivo y tener una mejor idea del aspecto de los datos que se están manejando, ya puede pasarse a la tarea de modelización teniendo una idea más precisa de lo que debe aplicarse. En este estudio fueron varios los métodos de inferencia que se aplicaron en los distintos campos de interés: el algoritmo de los k-vecinos más cercanos, árboles de decisión, bosques aleatorios, regresión lineal y regresión logística. Además, el ajuste de hiperparámetros se llevó cabo en aquellos métodos en los que era conveniente. 2.4.1. Algoritmo de los k-vecinos más cercanos También conocido como algoritmo KNN (de sus siglas en inglés: k-nearest neighbors) es un método de aprendizaje supervisado no paramétrico que se utiliza tanto en problemas de clasificación como de regresión. Las predicciones se basan en la similaridad del individuo 12 CAPÍTULO 2. PRELIMINARES 2.4. MÉTODOS DE INFERENCIA que se pretende predecir con otros individuos para los que el valor de la variable objetivo es conocida. El primer paso es seleccionar el número de vecinos que serán tenidos en cuenta (k) en base a la magnitud del conjunto con el que se esté trabajando y a su distribución. Generalmente este valor se obtiene comparando la eficiencia de distintos valores en una submuestra del conjunto total (conjunto de validación). Debe tenerse en cuenta también la distancia que va a utilizarse para medir la proximidad entre individuos. Generalmente la distancia utilizada para medir la similitud entre los individuos es la distancia euclídea, que equivale a la distancia de Minkowski cuando p=2, aunque cualquier exponente válido podría servir: d((x1, .., xm),(y1, .., ym)) =   m X j=1 |xj−yj|p  1/p Para cada uno de los individuos a predecir debe calcularse la distancia con todos los del conjunto de entrenamiento. Una vez se han medido todas las distancias, se seleccionan los k individuos de la muestra entrenamiento más cercanos al que se quiere predecir, son sus k vecinos más cercanos. En función de si se trata de un problema de clasificación o de regresión el procedimiento es distinto para el paso siguiente: Clasificación: se asigna la clase que presente mayoritariamente su conjunto de k vecinos más cercanos. Para evitar empates en el caso de tener solo dos categorías es conveniente escoger un valor impar para k. Regresión: se asigna el promedio de los valores de sus k vecinos más cercanos. Este proceso se repite para todos los individuos que se quieran predecir, no se construye un modelo realmente, por lo que la predicción de un nuevo caso puede suponer un coste computacional considerable para conjuntos de datos de gran tamaño. A cambio, se evita el coste asociado generalmente al proceso de aprendizaje, pero en muchos casos este coste puede ser asumible pero la aplicación del modelo requerirse en tiempo real por lo que se es más exigente con el tiempo de predicción. Más información en [19]. 2.4.2. Árboles de decisión Los árboles de decisión son una técnica de aprendizaje supervisado no paramétrico que puede utilizarse tanto para problemas de clasificación como de regresión. Es una técnica bastante popular ya que son sencillos de interpretar. Los árboles constan de una raíz desde la que se van generando de ramas y nodos. Los nodos internos presentan una división de los elementos del conjunto en función de las variables, y las ramas conectan los nodos entre sí. Si la variable utilizada para realizar la división es numérica, se realiza una división binaria; en caso de ser categórica puede hacerse también una división binaria o tener tantas alternativas como categorías haya en CAPÍTULO 2. PRELIMINARES 13 2.5. SOFTWARE registrados, y debía separarse en dos variables distintas. También pivot_wider fue útil para tratar con un conjunto en el que cuatro medidas distintas se registraban en una misma variable, habiendo una segunda variable para registrar de qué medida se trataba. Su complementario, pivot_longer, también se encuentra en este paquete. En caso de ser necesario también podría trabajarse con datos anidados gracias a nest yunnest, aunque no fue necesario para este proyecto en concreto. El paquete ggplot2 ([27]) es la herramienta por excelencia para crear gráficos a nivel profesional en este lenguaje. Trabaja con una estructura en capas, pudiendose añadir más y más detalles sobre un gráfico base, y pudiendo llegar a personalizar la estética con un enorme nivel de detalle de forma intuitiva gracias a su elegante gramática. Otras opciones para generar gráficos interactivos son plotly ohighchartr, entre otros, aunque ggplot2 sigue siendo el paquete más ampliamente utilizado para gráficos. El tratamiento de datos de tipo fecha fue sencillo gracias a lubridate, en este caso ocurre al igual que con readxl, aunque forme parte de tidyverse no se carga automáticamente al cargar el ecosistema . Este paquete está dedicado a poder leer fechas sin importar los distintos formatos que puedan presentar para así poder trabajar correctamente con ellas. Permite así extraer el mes o el año de una fecha sin tener que recurrir a expresiones regulares constantemente. Por último, el paquete forcats permite manejar las variables categóricas como factores, pudiendo reordenar sus niveles en función de otra variable (fct_reorder), en función de la frecuencia de cada categoría (fct_infreq), según una ordenación creada a mano (fct_relevel), o agrupar las cateogrías menos frecuentes en una única “otros” (fct_lump). 2.5.2. Ecosistema tidymodels y alternativas El ecosistema tidymodels (la información de este ecosistema y todos sus paquetes se encuentra en [15]) sigue la idea de tidyverse para combinar distintos paquetes y hacer que funcionen en conjunto, aunque en este caso la finalidad es simplificar y sistematizar la tarea de modelización y el aprendizaje automático. Este ecosistema fue creado por Hadley Wickham y Max Kuhn, entre otras personas, siendo el primero el creador de tidyverse, y el segundo el creador del paquete caret ([12]). Este ecosistema es precisamente la evolución natural de caret, y en varias ocasiones a lo largo de este trabajo se alterna entre ambas opciones como se describe más adelante. En Building predictive models in R using the caret package, [13], Max Kuhn habla sobre el uso de este paquete para construir modelos predictivos en R. También puede recurrirse a los vídeos de Julia Silge, o a su libro Tidy Modeling with R, [14], escrito junto con Max Kuhn para saber más del tema. Para comenzar, este ecosistema facilita la división en muestras de entenamiento y test con el paquete rsample y la función initial_split, de modo que ambas muestas queden almacenadas en un único objeto pero también puedan extrarse por separado con training ytesting. Esta función además permite indicar qué proporción de la muestra destinar al conjunto de entrenamiento (por defecto es 0.75) y también estratificar en función de una variable. De ese modo puede garantizarse que el valor de la variable objetivo sea similar en las muestras de entrenamiento y test. En caso de que no se realizara la división de forma estratificada y por algún casual los casos en ambos grupos tuvieran valores significativamente distintos, el ajuste se vería perjudicado pues los valores a predecir no 20 CAPÍTULO 2. PRELIMINARES 2.5. SOFTWARE son similares a aquellos utilizados para entrenar. Es por eso que es mejor estratificar de antemano, y al menos poder prevernir eso dentro del conjunto de datos disponibles. Sigue siendo posible, por supuesto, que en el futuro los datos que quieran predecirse sean muy distintos a aquellos con los que se entrenó y evaluó el modelo, lo cual no podría haberse previsto. En caso de empeorar mucho el rendimiento, tal vez obligase a generar un nuevo modelo para adaptarse a la nueva situación. Otro de los paquetes esenciales de tidymodels, el paquete parsnip, sirve para crear la base de los modelos, indicando de cuál se trata, si es de regresión o de clasificación, y añadiéndole un motor que se trata de otro paquete de R que ya está implementado (como puede ser rpart, [24]). El preprocesado de las variables se lleva a cabo con el paquete recipes, que incluye, entre otras muchas, funciones para: eliminar variables en caso de que tengan una correlación demasiado alta con step_corr. centrar y escalar las variables deseadas de forma simultánea, step_normalize, o bien solo centrar, step_center, o solo escalar, step_scale. imputar valores perdidos; hay una gran cantidad de funciones dedicadas a esto, pero en concreto la que fue usada para este trabajo fue step_impute_knn. eliminar automáticamente variables que tomen tomen solo un valor o tengan varianza muy cercana a cero mediante step_nzv. crear una categoría ficticia para las variables nominales, de modo que si se quiere predecir una observación que tenga un valor nunca visto en la muestra de entrenamiento no se produzca un error, step_novel. para las variables nominales con varias categorías, crear variables dummy que permitan traducir estas categorías en variables numéricas binarias, con step_dummy. Tanto los modelos como las recetas de preprocesado pueden combinarse en un flujo de trabajo con workflows, que puede después modificarse y actualizarse de manera sencilla. En caso de que se quiera mejorar el ajuste de hiperparámetros para el modelo, con el paquete tune puede realizarse esta tarea de forma eficiente. Finalmente, una vez se tenga establecido el modelo final, su eficiencia puede evaluarse con la libería yardstick, así como obtener las predicciones. Todo esto crea un flujo bastante bien definido para las tareas de modelización. Sin embargo, para el análisis que se llevó a cabo algunas de estas tareas fueron hechas de modo más manual debido a la familiaridad que se tenía con otras alternativas. En otro casos la eficacia de tidymodels fue lo suficientemente convincente como para optar por este ecosistema. Concretamente, gracias al paquete caret varias de estas tareas se pueden llevar a cabo sin un excesivo coste computacional. Para la división en dos muestras sí que se usó rsample aunque la función base sample también podría llevar a cabo esta división. Para los pasos posteriores, la función train de caret permite indicar el modelo a llevar a cabo con el argumento “method”, y aplicar un centrado y escalado de las variables numéricas con “preProcess”. También se puede hacer ajuste de hiperparámetros, aplicando por CAPÍTULO 2. PRELIMINARES 21 2.5. SOFTWARE ejemplo validación cruzada si se le indica, gracias a “tuneLength” si solo quiere indicarse el número de valores distintos para los parámetros que debe probar o “tuneGrid” si se quiere indicar expresamente los valores. Para los modelos mediante árboles de decisión y con bosques aleatorios se usaron respectivamente los paquetes rpart ([24]) o C50 ([20]) yrandomForest ([16]). 2.5.3. Otros paquetes Aparte de los ya mencionados, algunos paquetes no englobados en tidyverse fueron de gran utilidad. A la hora de generar el documento en cuestión se utilizó el paquete knitr ([33]), que es un motor dedicado a la generación de informes dinámicos con R. Gracias a esto se puede integrar código R en documentos Markdown y el uso de LaTeX para la estructuración. La estructura de este documento es la proporcionada por Luque-Calvo [18], este libro utiliza la propia plantilla para explicar la forma en que puede ser utilizada para llevar a cabo trabajos de fin de estudios con Rmarkdown. Para visualizar las tablas de forma más elegante se utilizó el paquete kableExtra ([34]), que está dedicado a simplificar la generación de tablas, permitiendo modificar su aspecto con sencillos argumentos. En el estudio de medidas de tendencia central y de dispersión de los datos fue útil el paquete skimr ([10]), así como para hacer un conteo de los casos perdidos en el caso de haberlos, ya que proporciona de forma automática tanto el número como el porcentaje de completitud. Para la creación de una aplicación interactiva que permite visualizar una gran cantidad de gráficas de forma sencilla y rápida se utilizó la librería shiny ([32]). Esta aplicación puede encontrarse en el apartado de contenido adicional, y gracias a ella es sencillo observar el comportamiento de algunas variables de este estudio sin tener que crear cada gráfica individualmente. Los manuales de los paquetes mencionados pueden encontrarse en la bibliografía en caso de querer consultarse de forma más extensas las funcionalidades de cada uno. 22 CAPÍTULO 2. PRELIMINARES Parte II Estudio estadístico de los datos 23 Capítulo 3 Importación, tratamiento y descripción de los datos Los datos utilizados para este estudio han sido recopilados por un equipo de ecólogos en el Parque Nacional de Doñana, Huelva. En dicha reserva hay una población considerable del animal de interés: el cangrejo rojo americano, Procambarus Clarkii. Este estudio se llevó a cabo entre 2004 y 2022, habiendo recopilado información sobre una gran cantidad de individuos de esta especie. Los datos del estudio, que fueron publicados el 11 de enero de 2023, constan de tres ficheros csv que a continuación serán analizados y ligeramente modificados para facilitar el tratamiento de la información. (Estos datos pueden consultarse en este enlace o en [5]). Además de eso, gracias a la colaboración entre estos ecólogos y el departamento de Ciencias de la Computación e Inteligencia Artificial de la Universidad de Sevilla en un proyecto en común, se dispone de otro conjunto de datos que no fue publicado. Estos datos han sido refinados en menor medida, y por tanto no se publicaron, pero al ser de una naturaleza distinta de aquellos que sí fueron publicados aportan una información muy valiosa, y por tanto se tomó la decisión de incluirlos también en este estudio. 3.1. Archivo OCC: ocurrencias individuales El archivo Don_ biom_red-swamp-crayfish _occ_20230111.csv contiene quince variables y 15483 registros. Sin embargo, la mayoría de las variables no aportan demasiada información al contener un mismo valor para todos los registros. Cada registro de este documento se refiere a una captura de un individuo, del cuál la información que se registra es: eventID: tiene el formato {ID localización}_{fecha AAAA-MM-DD}; esta información se separará en dos variables: IDL yfecha . occurrenceID: similar al anterior, pero añadiendo al final un contador que identifica de manera única a cada registro. Esto sirve para tener un identificador para cada cangrejo que fue estudiado, de aquí se extrae la variable Cangrejo. En el archivo “mof” se utiliza esta misma identifiación, por tanto pueden combinarse los dos conjuntos de datos. 25 3.1. ARCHIVO OCC: OCURRENCIAS INDIVIDUALES ocurrenceTime: tiene el formato H:M:S con un “z” al final de la cadena; se leerá ignorando esto último para anotar solo con el momento en que se registró la captura de cada individuo. individualCount: ya que cada registro corresponde a un individuo, siempre toma el valor 1 y por tanto no es una variable de interés para el estudio posterior. basicOfRecord: informa del método de captura; en este archivo solo toma el valor “Capture/Human Observation”; ya que todas se realizaron de la misma forma, tampoco tiene un interés para el estudio posterior. sex: toma los valores “Male” o “Female” en función de si el cangrejo era macho o hembra, aunque también hay una proporción considerable de datos perdidos (aproximadamente 16%) lifestage: toma los valores “Immature” o “Mature” en función de si se consideró que el cangrejo capturado era una cría o era ya adulto; en este caso, la cantidad de datos perdidos es mucho menor. Las siguientes variables se refieren a características de la especie, y al ser todos los registros sobre Procambarus Clarkii el valor que toman es el mismo en todo el conjunto; dicho valor se indica al lado de cada variable a continuación: •kingdom - “Animalia” •phylum - “Arthropoda” •class - “Malacostraca” •family - “Cambaridae” •genus - “Procambarus” •specificEpithet - “clarkii” •scientificName - “Procambarus clarkii” •scientificNameAutorshup - “Procambarus clarkii (Girard, 1852)” Por tanto, las variables que tienen relevancia para el análisis que se llevará a cabo son: ## IDL fecha Cangrejo occurenceTime ## Min. : 1.00 Min. :2004-06-09 Min. : 1 Length:15483 ## 1st Qu.: 17.00 1st Qu.:2008-04-22 1st Qu.: 3872 Class1:hms ## Median : 34.00 Median :2010-06-23 Median : 7742 Class2:difftime ## Mean : 54.37 Mean :2011-07-19 Mean : 7742 Mode :numeric ## 3rd Qu.: 55.00 3rd Qu.:2012-11-29 3rd Qu.:11612 ## Max. :490.00 Max. :2022-09-08 Max. :15483 ## sex lifestage ## Female:6504 Immature:8968 ## Male :6503 Mature :6504 ## NA’s :2476 NA’s : 11 ## ## ## 26 CAPÍTULO 3. IMPORTACIÓN, TRATAMIENTO Y DESCRIPCIÓN DE LOS DATOS 3.2. ARCHIVO MOF: INFORMACIÓN BIOMÉTRICA 3.2. Archivo MOF: información biométrica El archivo Don_ biom_red-swamp-crayfish _mof_20230111.csv contiene siete variables y 19599 registros. En este caso el formato es algo distinto al que se considera un conjunto de datos ordenados (cada fila un registro, cada columna una variable, cada celda un valor). Hay cuatro posibles medidas sobre las capturas individuales de cangrejos que en el archivo OCC se detallaban, pero para cada medida se ha generado una fila, en vez de tener una fila por cada cangrejo y cada medida en una variable. A continuación se detalla la información que contiene cada variable del archivo original, así como la posterior forma de organizarlo en el conjunto de datos que se utilizarán para el posterior análisis. occurrenceID: tiene el formato {ID localización}_ {fecha AAAA-MM-DD}_{ID cangrejo} por lo que lo que se separará en tres variables: IDL,fecha,Cangrejo. measurmentID: es igual que la anterior variable, pero añadiendo al final también un identificador que numera cada registro. Esto hace que sea un identificador único para el registro, pero de cara al análisis posterior no tiene utilidad, con las tres variables antes explicadas es suficiente. Nota: debería ser, al igual que en las siguientes variables, “measurement”, pero en el csv aparece con este nombre. measurementValue: variable numérica que indica la medida que se tomó, sin especificar unidades ya que eso se incluye en la siguiente variable. measurementUnit: para las variables de longitud el valor es “mm” pues se miden en milímetros y para la variable peso es “g” se mide en gramos. measurementType: toma cuatro valores distintos según la medida que se esté registrando: puede ser “total body length”, “cephalothorax length”, “cephalothorax width” o “weight” en función de si se está midiendo el largo del cuerpo, el largo del cefalotorax, el ancho del cefalotorax o el peso del cangrejo. measurementAccuracy: la precisión de la medida depende del tipo de medida: para el largo del cuerpo es 0.5 milímetros, para el largo y el ancho del cefalotorax es 0.1 milímetros, y para el peso es 0.2 gramos. measurementMethod: el método de medida depende del tipo de medida, toma solo un valor para cada una de las cuatro opciones y se especifica la forma en que se realizó la medición. Se puede ver que tres de las variables toman el mismo valor para cada uno de los cuatro posibles tipos de medida que se está registrando, por tanto es conveniente quedarse solo con el tipo de medida; de esta forma se reduce el número de variables sin perder información. Este conjunto será reordenado para que cada fila corresponda a un cangrejo. Se tendrán por tanto las variables fecha, IDL, Cangrejo, y las cuatro medidas (en muchos casos serán NA por no haberse tomado todas las mediciones en todos los cangrejos). Nota: Una vez realizada la ordenación, se puede observar que el cangrejo número 9, capturado el 4 de junio de 2010 en la localización con ID 9, tiene registrado el largo del cuerpo en 993 mm y ninguna otra medida. El siguiente valor más alto son dos cangrejos CAPÍTULO 3. IMPORTACIÓN, TRATAMIENTO Y DESCRIPCIÓN DE LOS DATOS 27 3.3. ARCHIVO EV: EVENTOS que miden 149 mm, por tanto se ha considerado que esta medición fue un error de medida o de registro, y se ha eliminado. El conjunto obtenido tras estas modificaciones tendría la siguiente forma: ## IDL fecha Cangrejo largo_cuerpo ## Min. : 1.00 Min. :2004-06-09 Min. : 1 Min. : 2.00 ## 1st Qu.: 17.00 1st Qu.:2008-04-22 1st Qu.: 3871 1st Qu.: 37.00 ## Median : 34.00 Median :2010-06-23 Median : 7742 Median : 63.00 ## Mean : 54.37 Mean :2011-07-19 Mean : 7742 Mean : 62.85 ## 3rd Qu.: 55.00 3rd Qu.:2012-11-29 3rd Qu.:11613 3rd Qu.: 86.00 ## Max. :490.00 Max. :2022-09-08 Max. :15483 Max. :149.00 ## ## largo_cefalotorax ancho_cefalotorax peso ## Min. : 1.30 Min. : 1.0 Min. : 1.00 ## 1st Qu.:11.70 1st Qu.:12.0 1st Qu.: 4.00 ## Median :18.70 Median :19.0 Median : 12.00 ## Mean :21.79 Mean :23.2 Mean : 16.93 ## 3rd Qu.:29.00 3rd Qu.:32.0 3rd Qu.: 26.00 ## Max. :78.00 Max. :68.0 Max. :117.00 ## NA’s :13552 NA’s :14434 NA’s :14344 Este conjunto de datos y el conjunto OCC pueden combinarse gracias al identificador de los cangrejos, de modo que en un único conjunto se tendrá la información sobre el sexo, estado de madurez, y las medidas de cada cangrejo, así como la fecha de captura, el lugar y su identificador. Como ambos conjuntos estaban bien relacionados, todos los registros que hay en uno los hay en el otro, así que a la hora de aplicar el join no hay que preocuparse sobre cuáles se mantienen, sino que todos están en el conjunto final llamado occ_mof. 3.3. Archivo EV: eventos El archivo Don_biom_red-swamp-crayfish_ev_20230111.csv contiene 18 variables y 444 registros, cada uno correspondiente a una localización geográfica y un día concreto, y se registran en ellos las capturas que se produjeron mediante trampas, no de forma individual como en los anteriores archivos. Hay muchas de estas variables que tienen un único valor, por tanto no se incluirán en el conjunto reducido que se usará en el análisis posterior. Concretamente, la información que se registra en las variables es: eventID: esta variable tiene el formato {ID localización}_{fecha en formato AAAA-MM-DD}, conviene que quede registrado en dos variables separadas, IDL yfecha, pero ya hay dos variables en este archivo que se corresponden con esa información, así que esta es redundante. intitutionCode: tiene un único valor: “ICTS-RBD”. Por tanto no aporta ninguna información relevante. 28 CAPÍTULO 3. IMPORTACIÓN, TRATAMIENTO Y DESCRIPCIÓN DE LOS DATOS 3.3. ARCHIVO EV: EVENTOS institutionID: también toma un único valor: “https://deims.org/bcbc866c-3f4f-47a8-bbbc-0a93df6de7b2”. datasetName: valor único: “ES_DONANA_Biom_red-swamp-crayfish_2004_2022_v020230101” eventDate: registra la fecha en formato AAAA-MM-DD, coincide con el valor que aparecía en la variable eventID. eventTime: registra la hora a la que se produjo la recogida de la trampa, el formato que tiene es H:M:S con un “z” al final, se selecciona únicamente la hora. continent: toma solo el valor “Europe”, las siguientes variables también tienen un único valor ya que todas las capturas de este archivo se corresponden con Doñana. country: solo toma el valor “ES”. stateProvince: solo toma el valor “Huelva”. location: solo toma el valor “Doñana”. localityID: indica el identificador numérico de la localización en la que se produjeron las capturas, se corresponde con la primera parte de la variable eventID. locality: para cada localización aparte de un código numérico se asigna un nombre para identificarla. Es una información útil, pero para el análisis posterior esta variable no aporta nada que no se conozca con el identificador numérico, así que se eliminará. decimalLatitude: coordenada Latitud de la localización en la que se produjeron las capturas medida en grados decimales WGS 84. decimalLongitude: coordenada Longitud de la localización de las capturas medida en grados decimales WGS 84. habitat: indica si la localización se trata de una laguna o una marisma. sampleSizeValue: variable numérica que indica el número de capturas que se registraron. sampleSizeUnit: solo toma el valor “trap”, todas las capturas se realizan de la misma forma. sampleSizeEffort: solo toma el valor “24 hours”. Una vez se eliminan las variables que toman un único valor o no aportan nada adicional, el conjunto anterior se reduce considerablemente. También se eliminan aquellos registros que tengan perdido el valor de SampleSizeValue, ya que no aportan nada. El aspecto final del conjunto de datos es el siguiente: CAPÍTULO 3. IMPORTACIÓN, TRATAMIENTO Y DESCRIPCIÓN DE LOS DATOS 29 4.1. FASE EXPLORATORIA Figura 4.1: Captura documento shiny La segunda pestaña, Variables ambientales, tiene un selector que filtra los datos del conjunto pres_abs para reducirse a los de una localización en concreta. Hay también un selector para cada variable ambiental, cada uno añade su correspondiente variable a la gráfica en caso de ser seleccionado, y la elimina en caso de deseleccionarse. La tercera pestaña, Proporción de individuos, permite ver información sobre los cangrejos capturados en una localización concreta y para un rango de fecha personalilzable con un selector. Hay dos subpestañas, en función de si quiere verse la proporción de Machos/Hembras o de individuos Maduros/Inmaduros. Se muestra el número de capturas por cada fecha con algún registro, así como la proporción total y en cada fecha de la variable de interés. La cuarta pestaña, Tamaño cangrejos, muestra un histograma para cada una de las cuatro variables relacionadas con el tamaño que se tomaron en los cangrejos. El número de barras se puede elegir con un selector, así como la localización y el rango de fechas sobre los que quiere visualizarse la información. La quinta pestaña, Modelo reactivo, sirve para visualizar los resultados de los modelos KNN de clasificación en función de determinados parámetros. Los selectores de esta pestaña son: la variable objetivo (Sexo o Madurez), el número de vecinos utilizado y el exponente de la distancia de Minkowski utilizado. El modelo se genera al pulsar el botón de “Actualizar”, y se muestran la precisión, el área bajo la curva roc, y una gráfica mostrando la matriz de clasificación. Como se comentó brevemente en el apartado de ajuste de hiperparámetros, la rapidez que permite esta aplicación para variar los parámetros y ver cómo eso afecta a los resultados sirve para resaltar la importancia de un buen ajuste. La aplicación permite generar un gran número de gráficas disintas en un tiempo reducido, pero estas no siempre aportan información relevante sobre el conjunto de datos. Es por eso que algunos de los puntos más destacables que se descubrieron a lo largo de la fase exploratoria se detallan en los apartados siguientes de este capítulo. Puede verse en la Figura 4.1 una captura de la aplicación generada. 36 CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 4.2. DIFERENCIAS ENTRE HÁBITATS 4.2. Diferencias entre hábitats Se pueden representar las coordenadas de las distintas localizaciones, obteniendo así una visión orientativa a modo de mapa del terreno sobre el que se realizó el estudio. Incluso puede identificarse para cada localización si se trata de una Laguna o de una Marisma. 36.9 37.0 37.1 −6.7 −6.6 −6.5 −6.4 −6.3 Longitud decimal Latitud decimal habitat Lagunas temporales Marismas Figura 4.2: Localizaciones geográficas Una vez visto que ambos tipo de hábitat están en zonas bastante diferenciadas, podría surgir la duda de si en cada una de estas zonas el tamaño de los cangrejos tiene una diferencia significativa o si esta diferencia en la localización no es un factor relevante. 0 50 100 150 Lagunas temporales Marisma Habitat Largo Cuerpo Figura 4.3: Longitud de los cangrejos según el hábitat CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 37 4.3. VARIABLES AMBIENTALES A simple vista parece que la distribución que sigue la variable “Largo cuerpo” en ambos hábitat es bastante similar, los cuartiles tienen valores muy similares en ambos grupos. Esto permitirá que los modelos que utilicen las variables relacionadas con el tamaño no tengan que ser tratados de forma distinta, sino que pueda englobarse la población total. 4.3. Variables ambientales Es importante observar las correlaciones entre las variables, ya que de cara a algunos modelos podría suponer un problema en caso de haber multicolinealidad. Dado que en el conjunto de datos ambientales hay quince variables numéricas, se podría representar su matriz de correlaciones, pero el número de pares posibles es excesivamente elevado, y la mayoría de valores son tan bajos que no merece la pena. En vez de hacer la representación total, se han extraído los dos pares de variables que tienen una correlación en valor absoluto mayor de 0.8: Salinidad-Conductividad y Oxígeno disuelto en porcentaje-Oxígeno disuelto en valor absoluto. Para ambos pares se hace la representación de los pares de valores colocando cada variable en un eje, así como el ajuste de la variable en el eje Y respecto a la del eje X mediante geom_smooth(). 0 10 20 30 40 0 10 20 Salinidad Conductividad Conductividad vs Salinidad 0 10 20 0 50 100 150 200 250 Oxígeno disuelto en porcentaje Oxígeno disuelto en valor absoluto Oxígeno en valor absoluto vs Oxígeno en porcentaje Figura 4.4: Variables ambientales más correlacionadas La relación que se observa entre el primer par de variables se debe a propiedades químicas que eran esperables, la salinidad del agua favorece su conductividad. Viendo el ajuste se puede observar que ambas siguen una tendencia casi lineal con una dispersión considerablemente pequeña, apenas hay valores que se alejen de la diagonal. La relación entre el segundo par de variables se debe a que ambas tratan sobre el mismo tipo de dato, solo que en distintas escalas. Incluso así, el ajuste que se puede realizar entre ambas no es tan bueno como en el primer caso, conforme aumenta el valor de las variables también aumenta la dispersión. La correlación que tienen ambos pares de variables es muy alto, lo cuál significa que de cara a futuros modelos se debe tener en cuenta que incluir ambas puede producir problemas derivados de la multicolinealidad. Concretamente, del primer par se eliminará la Salinidad, ya que tiene 237 valores perdidos, mientras que la Conductividad tiene 11. De ese modo se evita inutilizar un gran número de registros, pero sin perder apenas información. 38 CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 4.3. VARIABLES AMBIENTALES Aparte de estos dos pares de variables, las únicas con una correlación considerable eran Fecha-Año, ya que Fecha fue obtenida de la combinación entre Año y Mes, pero su finalidad es únicamente para que en las gráficas pueda apreciarse mejor la evolución de las variables. Por otro lado, la Temperatura está algo relacionada (aunque bastante menos que en los casos anteriores) con la Fecha y con el Año, en cambio no con el Mes. La relación que tienen estas variables no es tan clara como las dos anteriores, pero también es conveniente ver cómo evoluciona. La temperatura tiene un comportamiento que parece ser cíclico a lo largo de los años del estudio; si se ignora el año y se considera solo la variable Mes, puede verse que esta evolución es la que cabía esperar: mayores valores en los meses centrales (verano) y valores más reducidos en el invierno. 0 10 20 30 40 2004 2007 2010 2013 2016 Fecha Temperatura Temperatura según Fecha 0 10 20 30 40 1 2 3 4 5 6 7 8 9 10 11 12 Mes Temperatura Figura 4.5: Evolución de la temperatura Por otro lado, es relevante observar cómo varían los valores de Conductividad (y por extensión los de Salinidad) en función de si se trata de una Marisma o una Laguna. 0 10 20 30 40 2008 2010 2012 2014 Fecha Conductividad Habitat Lagunas temporales Marisma Figura 4.6: Conductividad según hábitat Puede observarse que en las Marismas se dan los valores más altos de esta variable, aunque también se dan algunos valores reducidos, mientras que en las Lagunas los valores CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 39 4.4. TAMAÑO SEGÚN SEXO Y MADUREZ son bajos en la mayoría de los casos. Esto indica que, de cara a los análisis que engloben esta información en distintas localizaciones, se debe hacer una separación entre las Lagunas y las Marismas. 4.4. Tamaño según sexo y madurez Es interesante plantearse también si el tamaño de los cangrejos varía mucho en función de si son macho o hembra. Esta información podría ser útil para poder clasificar el sexo de los cangrejos en función de su tamaño. largo_cuerpo peso ancho_cefalotorax largo_cefalotorax Female Male NA Female Male NA Female Male NA Female Male NA 0 20 40 60 80 0 30 60 90 120 0 20 40 60 0 50 100 150 Sexo Figura 4.7: Distribución de las medidas de tamaño según sexo del cangrejo La distribución entre machos y hembras parece ser similar en las cuatro variables, pero en cambio los valores perdidos actúan de manera bastante diferenciada. Esto puede dificultar la clasificación de los valores perdidos en una de las dos categorías si se intentase realizar únicamente en base a estas cuatro variables. Para entender el motivo de que la distribución de los valores perdidos tenga este comportamiento, puede ser útil realizar esta misma representación pero respecto a la madurez de los cangrejos. 40 CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 4.5. CAPTURAS EN TRAMPAS Tabla 4.1: Tabla cruzada Madurez y Sexo Female Male NA Immature 3438 3079 2451 Mature 3059 3419 25 NA 6 5 0 largo_cuerpo peso ancho_cefalotorax largo_cefalotorax Immature Mature NA Immature Mature NA Immature Mature NA Immature Mature NA 0 20 40 60 80 0 30 60 90 120 0 20 40 60 0 50 100 150 Sexo Figura 4.8: Distribución de las medidas de tamaño según madurez del cangrejo En este caso las variables sí tienen comportamientos significativamente distintos entre los individuos Maduros e Inmaduros, lo cuál es bastante razonable. Significa también que una clasificación sobre la etapa de Madurez del individuo probablemente fuese más adecuada que una sobre el sexo del individuo. La clave está en observar qué relación siguen los valores perdidos de cada una de las variables como puede verse en la Tabla 4.1. Puede apreciarse que casi todos los datos perdidos del sexo son en individuos inmaduros, seguramente porque esto dificulta su clasificación, y al tener un tamaño menor que los individuos maduros eso hace que los valores perdidos del sexo tengan un tamaño considerablemente menor que los que fueron catalogados como machos o como hembras. 4.5. Capturas en trampas En un principio esta variable aparentaba tener una gran cantidad de información disponible, lo cual es cierto, pero la forma en que fueron recogidos dificulta la interpretación del tamaño de la población de cangrejos en base a estos datos. No se trata de una variable que sea constante en el tiempo, sino que en cada fecha el número de localizaciones en que se CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 41 4.5. CAPTURAS EN TRAMPAS realizaron capturas es disinto. Tampoco los registros son periódicos, sino que se reparten a lo largo de los años del estudio sin seguir un patrón y siendo escasos en algunos períodos extendidos de tiempo. El número de registros en una localización concreta tampoco son suficientemente numerosos como para poder estudiar la evolución de los cangrejos en ninguna de las localizaciones disponibles. Todo esto hace que interpretar esta variable para sacar conclusiones sobre el tamaño de la población en Doñana a lo largo del estudio sea complicado. Por tanto, de cara a tratar de estimar la evolución que ha seguido la población de cangrejos, es difícil representarlo gráficamente. Tal vez un mes hubo menos capturas, pero porque solo se visitaron un pequeño número de localizaciones. En cambio, si en otro mes que se visitaron muchas más localizaciones o a lo largo de más días se hubiera alcanzado el mismo número de capturas. Es necesario poder ver reflejado que, en esa situación, en cada lugar se habrían recogido menos cangrejos (aunque de forma acumulada formasen un número considerable). La forma de corregir esto será mediante una transformación de los datos, pero antes de representarlo conviene mostrar los valores obtenidos sin importar fecha ni localización, para tener una idea del rango que abarcan y la frecuencia de aparición. 0 50 100 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Capturas en trampas Frecuencia Figura 4.9: Distribución capturas en trampas Parece ser que la distribución seguida es bimodal, siendo los dos valores con mayores frecuencias 5 y 9 individuos por trampa y día. Cabe destacar que la gráfica anterior engloba los dos tipos de hábitat, pero cada uno de ellos tienen una distribución considerablemente distinta. 42 CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 4.5. CAPTURAS EN TRAMPAS Lagunas temporales Marismas 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 0 30 60 90 120 Capturas en trampas Frecuencia Figura 4.10: Distribución capturas en trampas según hábitat En las lagunas temporales la frecuencia más habitual es 5, mientras que en las marismas es 9 (aunque siga habiendo ligeros picos en los valores contrarios). Esto indica que en cada tipo de hábitat, o tal vez incluso en cada localización distinta, los valores de cangrejos capturados pueden ser considerablemente distintos. Por tanto, se confirma que el detalle antes mencionado de separar los datos según el hábitat debe ser tenido en cuenta si se quiere analizar la evolución de las capturas a lo largo del tiempo. Para poder hacer este análisis y tener en cuenta al mismo tiempo la localización en que se registraron las capturas, se realizará una poderación de los valores antes mostrados. La forma de hacerlo será evaluando la media de capturas en cada una de las localizaciones, y respecto a esa media se ponderarán las capturas obtenidas cada día por cada localización. De ese modo, un valor cercano a 1 representaría un número de capturas estándar, cuanto mayor fuese el valor este mayor sería el número de capturas en ese día respecto a lo que se suele encontrar, y viceversa. CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 43 4.5. CAPTURAS EN TRAMPAS 0.5 1.0 1.5 2.0 2.5 2005 2010 2015 2020 Fecha Capturas ponderadas Figura 4.11: Capturas ponderadas Puede verse que en general los puntos están distribuidos de forma bastante equilibrada. En las fechas que en algunas localizaciones se realizaron más capturas hubo otras en las que se realizaron menos, lo cual parece indicar que se trataba de una fluctuación local, pero no de una evolución por parte de la población local de cangrejos. Por ver de forma más concreta el comportamiento que seguían las capturas, pueden visualizarse los seis localizaciones que tuvieron mayor número de capturas: 0.5 1.0 1.5 2.0 2.5 2005 2010 2015 2020 Fecha Capturas ponderadas Localización 9 15 28 30 31 34 Capturas ponderadas en localizaciones de mayor estudio 4 8 12 2005 2010 2015 2020 Fecha Capturas por trampa habitat Lagunas temporales Marismas Localización 9 15 28 30 31 34 Capturas por trampa en localizaciones de mayor estudio Figura 4.12: Capturas en localizaciones más estudiadas En la gráfica de la izquierda se pueden ver las capturas ponderadas, todas están alrededor del valor 1 aunque pueden observarse fluctuaciones a lo largo del tiempo. En la gráfica de la derecha se representan los valores reales, además de diferenciar con el tipo de línea aquellas localizaciones que son Marismas de las que son Lagunas. Puede verse que los valores se agrupan sobre todo a alturas 5 y 9, que eran los que se vio que eran más comunes, aunque a diferencia de lo que cabría esperar los dos hábitats se encuentran mezclados, no hay una división clara. Una de las marismas tiene valores más cercanos a 44 CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 4.5. CAPTURAS EN TRAMPAS 5, que era habitual en las lagunas, y una de las lagunas tiene valores más cercanos a 9 la mayoría del tiempo. Algo que llama la atención de ambas gráficas es que en la localización 34 los dos últimos registros reflejan una disminución drástica de la cantidad de cangrejos capturados, a pesar de que en los anteriores registros casi siempre tenía valores altos. Esto podría tener múltiples explicaciones, ya sea por un desplazamiento de los cangrejos de esta localización a alguna cercana, o debido a las condiciones de dicho lugar. Sin embargo responder a esta cuestión quedaría en mano de aquellos ecólogos que realizaron el trabajo de campo. Otra pregunta que puede hacerse es si la cantidad de cangrejos se ve influenciada por la estación. Para tratar de responder a la misma, se agrupan las capturas en intervalos de tres meses y se suman todas las que se realizaron por cada localización. Después de eso vuelve a aplicarse la ponderación, pero con los nuevos valores de capturas. 0.5 1.0 1.5 2.0 2.5 2005 2010 2015 2020 Fecha agrupada por estación Capturas ponderadas Figura 4.13: Capturas ponderadas agrupadas por estación La distribución agrupando las capturas por trimestre no varía demasiado, pero al reducir el número total de puntos en la gráfica ayuda a que se vea más claramente y puede representarse la media de cada fecha sin que sea demasiado confisa la gráfica. De nuevo, no parece que haya ido evolucionando por épocas, sino que en cada estación hubo valores por encima y por debajo de la media habitual en distintas localizaciones. Los mayores picos se corresponden a fechas con pocos registros (uno o dos en los casos más altos) aunque el segundo trimestre en el que se realizó el muestreo sí que presenta sus ocho registros con valores por debajo de la media. Lo mismo ocurre al visualizar las seis localizaciones con mayor número de registros. En algunas las capturas se mantienen de forma más estable mientras que en otras está más disperso, pero en general el comportamiento no dista demasiado del que se vio en la Figura 4.12 con las capturas diarias. El IDL 28 es el único que muestra algunos valores elevados que no tenía con la anterior agrupación y ponderación. CAPÍTULO 4. ANÁLISIS DESCRIPTIVO 45 Capítulo 5 Introducción a los modelos predictivos Tras haber realizado el estudio descriptivo de los datos, ya se tiene una idea general de cómo se comportan y cómo interaccionan entre ellas las variables, por lo que ha llegado el momento de intentar crear modelos que las ajusten de la mejor forma posible. Los modelos predictivos pueden ser de dos tipos: de clasificación, si tratan de predecir una variable categórica, o de regresión, si tienen como objetivo la predicción de una variable numérica. En este capítulo se hará un esbozo del guion que se seguirá para crear los modelos, con algunos trozos del código para explicar cómo fueron llevados a cabo. Los otros dos capítulos que componen este bloque se dedican a explicar dichos modelos, reflejando los resultados que se obtuvieron para cada una de las variables, y las conclusiones que se alcanzan. En cada uno de los capítulos se estudian varias variables, aplicando distintos métodos a algunas de ellas para poder compararlos entre ellos, y también hacer comparaciones entre las variables. En el de clasificación son la variable Sexo y la variable Madurez las que se intentan predecir, utilizando como variables predictoras las variables de tamaño (largo del cuerpo, largo del cefalotórax, ancho del cefalotórax y peso), aunque en ocasiones se usan todas y en otras solo el largo del cuerpo. En el de regresión se estudian las variables de tamaño, intentando predecirlas en función de las demás, y la variable de las capturas ponderadas, en función de las variables ambientales. 5.1. Ideas iniciales Antes de explicar los modelos llevados a cabo, merece la pena mencionar aquellas ideas que no llegaron a incorporarse en este estudio, y el motivo de que no fueran incluidas en el producto final. En aquellas variables que trataban sobre individuos concretos (tamaño, sexo, madurez) la línea a seguir estaba bastante clara, lo lógico era intentar modelizar cada una de ellas en función del resto de las variables. No es coherente intentar predecir lo que medirá un hipotético cangrejo futuro aleatorio, sin tener ni idea de alguna otra variable de tamaño. Tampoco clasificar si será maduro o inmaduro sin saber nada de lo que mide. Por tanto, aunque como se verá más adelante los modelos no siempre sean tan buenos como se pretendía, con los datos disponibles era lo único que podía hacerse, o al 53 5.1. IDEAS INICIALES menos lo único que humildemente pudimos concebir. La cuestión a ir desarrollando fue la elección de los modelos, y el ajuste de hiperparámetros en cada uno de ellos. En cuanto a las variables ambientales, el asunto se complicaba de forma considerable. Por la forma en que se toman las mediciones (de forma mensual) una primera idea podría ser tratar esos datos como series temporales, ver su evolución a lo largo del tiempo e intentar modelizar eso para alguna variable como la temperatura o el pH. Sin embargo, la forma en que están tomados los datos hizo que eso fuera impracticable. Las mediciones se tomaban en localizaciones distintas, en los meses que había al menos un registro, generalmente había registros sobre varias localizaciones distintas. Pero en un gran porcentaje de los meses no se visitaba ninguna de las localizaciones. Intentar estudiar una localización en concreto era imposible, pues los registros estaban muy espaciados en el tiempo, y la cantidad de valores perdidos en medio era superior a la cantidad de registros tomados. Las técnicas habituales para tratar con series temporales requieren que no haya valores perdidos, y en caso de haberlos lo apropiado sería intentar imputarlos, pero no es coherente imputar un número de datos mayor que la cantidad disponible, sobre todo siendo esta tan pequeña. Por tanto, predecir localizaciones individuales quedó descartado, así que habría que agruparlas de alguna forma para intentar tener suficientes registros como para crear un buen modelo. Intentar agrupar aquellas que fueron medidas en una misma fecha tenía como inconveniente que la mayoría de variables tenían un amplio rango y generalmente en un mismo mes podían alcanzar un rango prácticamente igual al rango total de la variable. Para verlo de forma más concluyente, se calcularon los valores máximos y mínimos alcanzados por cada variable en una fecha, además el rango de dicha fecha y la desviación típica de las observaciones. Se muestran dichos datos para las variables Temperatura y Nitrato, que tienen escalas bastante diferentes y en ambos hábitats son similares, pero en ninguno de los casos parece que agrupar los datos medidos en el mismo mes vaya a ser adecuado. pres_abs_Temp <- pres_abs %> % filter(!is.na(Temperatura)) %> % group_by(Fecha) %> % summarise( across( c(Temperatura), list(max = ~max(., na.rm = TRUE), min = ~min(., na.rm=TRUE), rango = ~(max(., na.rm = TRUE)-min(., na.rm = TRUE)), sd = ~sd(., na.rm=TRUE)), .names = "{.col}_{.fn}" ) ) %> % ungroup() %> % select(-Fecha) summary(pres_abs_Temp) ## Temperatura_max Temperatura_min Temperatura_rango Temperatura_sd ## Min. : 5.63 Min. : 0.61 Min. : 0.000 Min. :0.1131 ## 1st Qu.:15.06 1st Qu.: 8.40 1st Qu.: 3.250 1st Qu.:1.6683 ## Median :19.46 Median :13.17 Median : 6.630 Median :2.2641 54 CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 5.1. IDEAS INICIALES ## Mean :19.23 Mean :13.26 Mean : 5.972 Mean :2.2612 ## 3rd Qu.:23.91 3rd Qu.:18.53 3rd Qu.: 8.140 3rd Qu.:2.8410 ## Max. :38.90 Max. :25.95 Max. :15.900 Max. :5.0727 ## NA’s :6 pres_abs_nit <- pres_abs %> % filter(!is.na(Nitrato)) %> % group_by(Fecha) %> % summarise( across( c(Nitrato), list(max = ~max(., na.rm = TRUE), min = ~min(., na.rm=TRUE), rango = ~(max(., na.rm = TRUE)-min(., na.rm = TRUE)), sd = ~sd(., na.rm=TRUE)), .names = "{.col}_{.fn}" ) ) %> % ungroup() %> % select(-Fecha) summary(pres_abs_nit) ## Nitrato_max Nitrato_min Nitrato_rango Nitrato_sd ## Min. : 0.104 Min. : 0.048 Min. : 0.000 Min. : 0.1773 ## 1st Qu.: 5.204 1st Qu.: 1.102 1st Qu.: 3.483 1st Qu.: 2.8190 ## Median : 67.688 Median : 1.994 Median : 40.128 Median : 17.6047 ## Mean :177.884 Mean : 23.641 Mean :154.243 Mean : 60.7002 ## 3rd Qu.:188.350 3rd Qu.: 9.052 3rd Qu.:149.202 3rd Qu.: 76.3773 ## Max. :998.500 Max. :463.700 Max. :993.293 Max. :324.3659 ## NA’s :6 Puede verse en el rango de la temperatura que en un mismo mes puede haber una diferencia de casi 16 grados en las mediciones de dos localizaciones distintas, por lo que intentar agruparlos con el valor medio daría lugar a un error de unos 8 grados, lo cual es demasiado. Para el nitrato, el rango máximo alcanzado es de 993, cuando el rango total de la variable no llega a 1000, por lo que en un mismo mes se han registrado valores cercanos a lo más bajo y a lo más alto que llega a alcanzar la variable en todo el estudio. Se descarta también la posibilidad de utilizar el valor medio de un mes para así tener valores en un número de meses mayor. Por tanto, con los datos disponibles no se puede llevar a cabo un estudio con el enfoque de series temporales. Se probaron también técnicas de data augmentation como muestreo con reemplazamiento para así tener conjuntos de mayor tamaño, pero no se obtuvo ninguna mejora en los modelos. Finalmente, lo que pudo remediarse fueron los datos perdidos en algunos registros, aplicando imputación mediante knn, ya que de esa forma aquellos registros a los que les faltasen una o dos variables no tenían que ser descartados automáticamente. Esto se hizo mediante la función step_impute_knn de tidymodels tal como se puede ver en la sección 5.4. CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 55 5.2. ANÁLISIS DISTRIBUCIONAL También se hizo el intento de utilizar redes neuronales para intentar crear un modelo que mejorase los que se habían creado. Se partió de una red sencilla, utilizando la librería neuralnet, para realizar la clasificación de la variable Sexo. Tanto con la variable “Largo cuerpo” como con las cuatro variables del tamaño se obtenían resultados similares a los que se verán en la sección 6.1. Estos modelos están en auge dada su versatilidad y capacidad de predicción, pero es cierto que requieren de un gran entrenamiento y una dedicación considerable para obtener un modelo que realmente aproveche todo lo que esta técnica pueda ofrecer. Finalmente, terminó descartándose abarcar ese ámbito en este trabajo, ya que no era trivial conseguir un modelo que mejorase los que ya se habían desarrollado. 5.2. Análisis distribucional Según el manual sobre Tendencias e índices de seguimiento de datos, conocido como TRIM (TRends & Indices for Monitoring data) [21], dedicado precisamente a analizar datos de esta naturaleza, el conteo de poblaciones debería seguir una distribución Poisson. Por tanto es lo que utilizan para intentar estimar poblaciones en dicho manual. Sin embargo, también inciden en que hay que tener en cuenta que en los casos reales la varianza suele ser mayor a lo que se esperaría de una Poisson. Aunque el software del que hablan en su manual, TRIM, no es el que se utilizará para este estudio, sigue siendo de interés comprobar si la información sobre el conteo de los cangrejos puede aproximarse con esta distribución. En las gráficas antes representadas pudo apreciarse que la distribución que siguen las capturas con trampas en cada hábitat tenían una distribución algo distinta. Para saber si esta diferencia es realmente significativa, se realizará el análisis de la varianza ANOVA. ## Df Sum Sq Mean Sq F value Pr(>F) ## ev$habitat 1 328 328.0 63.72 1.63e-14 *** ## Residuals 389 2002 5.1 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 La hipótesis de igualdad de medias de la variable entre ambos conjuntos es rechazada (p-valor del orden de 10−14), por tanto se tratará con la información de ambos hábitats por separado. De modo que el análisis sobre bondad de ajuste de la distribución Poisson se realizará por separado para ambos subconjuntos. ## ## Goodness-of-fit test for poisson distribution ## ## X^2 df P(> X^2) ## Pearson 275.1126 13 3.802528e-51 En el caso de las marismas se rechaza la hipótesis de que los datos sigan una distribución Poisson (p-valor del orden de 10−51), lo cual significa que un ajuste realizado con el software TRIM no sería apropiado, y tampoco permite hacer estimaciones con este tipo de ajuste. 56 CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 5.3. MODELOS DE CLASIFICACIÓN ## ## Goodness-of-fit test for poisson distribution ## ## X^2 df P(> X^2) ## Pearson 271.5248 14 9.953403e-50 También para el caso de las lagunas temporales la hipótesis nula es rechazada (p-valor del orden de 10−50) por lo que se concluye que los datos no siguen una distribución Poisson. Ya que no cumplen lo que se esperaba de una variable de conteo poblacional, se puede intentar estimar con una distribución normal. Sin embargo, esto tampoco tiene éxito, ya que aplicando ambos conjuntos el test de normalidad de Shapiro-Wilk por separado se obtienen p-valores del orden de 10−12 (marismas) y 10−15 (lagunas). ## ## Shapiro-Wilk normality test ## ## data: . ## W = 0.86173, p-value = 4.213e-12 ## ## Shapiro-Wilk normality test ## ## data: . ## W = 0.81464, p-value = 9.156e-15 Por tanto se rechaza que la distribución del número de capturas de cangrejos realizadas en trampas siga una distribución normal ni en las Marismas ni en las Lagunas. El enfoque paramétrico no es apropiado para los datos disponibles, por lo que se intentará recurrir a las variables ambientales para ver si esa información adicional ayuda a garantizar un mejor ajuste de la variable. 5.3. Modelos de clasificación Para los modelos de clasificación se trató con las variables Sexo y Madurez. La primera de ellas tiene un porcentaje mucho mayor de valores perdidos, lo cual supone que la muestra disponible es menor. En ambos casos se utilizó un 70% de los registros completos para la muestra de entrenamiento y un 30 % para la muestra test, de forma estratificada, mediante el comando initial_split: occ_mof_comp_sexo_split <- occ_mof_completo_sexo %> % initial_split(prop=0.7,strata = sex) occ_mof_comp_sexo_ent <- occ_mof_comp_sexo_split %> % training() occ_mof_comp_sexo_test <- occ_mof_comp_sexo_split %> % testing() CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 57 5.3. MODELOS DE CLASIFICACIÓN Se aplicaron tres métodos de clasificación a cada una de las variables: k-vecinos más cercanos, bosques aleatorios y regresión logística. Para el caso de los k-vecinos más cercanos, los modelos se construyeron mediante la librería caret, y se realizó el ajuste del parámetro k mediante validación cruzada. El número de valores comprobados fueron 11, como puede verse en el código, y antes de realizar el modelo se aplicó un centrado y escalado de las variables predictoras, de modo que todas las variables mantuvieran a priori un peso similar a la hora de realizar el cálculo de distancias. ctrl <- trainControl(method="cv",classProbs=TRUE, summaryFunction = defaultSummary ) KNN_todas_sexo <- train(sex ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_sexo_ent, method = "knn", trControl = ctrl, preProcess =c("center","scale"), tuneLength=11 ) Este mismo procedimiento se repitió utilizando como variable predictora únicamente el largo del cuerpo, y lo mismo para modelar la variable Madurez. Como en bosques aleatorios no tiene sentido plantearse el modelo solo según el largo del cuerpo, esto se reservó para crear un árbol de decisión mediante la librería C50: arbol_largo_sexo <- C5.0(sex ~ largo_cuerpo, data = occ_mof_largo_sexo_ent) El bosque aleatorio se hizo con la librería randomForest. El parámetro mtry se fijó en 3, sin aplicar ajuste, ya que se consideró que la mejora que se obtenía no era suficiente como para justificar el coste computacional. bosque_todas_sexo <- randomForest(sex ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_sexo_ent, mportance = TRUE, replace = TRUE, ntree=100,mtry=3) La regresión logística se realizó mediante el comando glm, en este caso no fue necesario aplicar ingeniería de características. logi_todas_sexo <- glm(sex ~ largo_cuerpo + ancho_cefalotorax + largo_cefalotorax + peso, family = binomial, data = occ_mof_comp_sexo_ent) 58 CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 5.4. MODELOS DE REGRESIÓN El código reflejado en esta sección tiene como objetivo dejar claros los procedimientos seguidos, para así en el capítulo siguiente poder mostrar únicamente los resultados de forma más resumida. El código utilizado en cada uno de los modelos se encuentra en el Apéndice A. 5.4. Modelos de regresión Las variables que fueron modelizadas mediante modelos de regresión fueron el largo de cefalotórax, ancho del cefalotórax, peso y las capturas ponderadas. En este caso, el porcentaje dedicado a la muestra de entrenamiento fue el 80%, pues para la última variable el número de registros es reducido y era necesario un porcentaje algo superior para entrenar los modelos. La división se generó de nuevo con initial_split de forma estatificada. Por cuestiones que se explican de forma más extensa en el capítulo de Modelos de regresión, para el ancho y el largo de cefalotórax se realiza también una división del conjunto inicial en función de la fecha de captura de los cangrejos. Una vez creada esa división, se crearon modelos de regresión lineal mediante lm para cada subconjunto y cada variable. El error se midió mediante validación cruzada, aplicando 10 pliegues, mediante la función cv.lm tal como se ve en el código a continuación: reg_mof_larcf_f1 =lm(largo_cefalotorax ~ largo_cuerpo, data=training(mof_larcf_f1_split)) VC_mof_larcf_f1 <- cv.lm(largo_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_larcf_f1_split), seed = 1909,plotit = FALSE, printit = FALSE) Para la variable peso los resultados obtenidos mediante regresión lineal no eran suficientemente buenos, ni realizando la división que se hizo para las otras dos variables. Por la forma de las predicciones del modelo, se realizó la transformación logarítmica, y gracias a eso mejoraron los resultados de forma notable. De nuevo el error cometido se midió mediante validación cruzada. reg_log_peso =lm(log(peso) ~ largo_cuerpo, data=training(mof_pes_split)) Como esta variable no había obtenido buenos resultados para la regresión lineal sin aplicar la transformación, se realizó también el ajuste mediante bosque aleatorio utilizando las tres otras variables de tamaño como predictoras. En este caso se llevó a cabo con la librería caret, y se aplicó el ajuste del parámetro mtry pudiendo tomar valores 1, 2 o 3. Con esto pudo obtenerse también un buen ajuste. bosque_peso_tuning <- train(peso ~ ancho_cefalotorax + largo_cefalotorax + largo_cuerpo, data = training(mof_comp_split), method = 'rf',tuneGrid = expand.grid(.mtry = 1:3), CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 59 5.4. MODELOS DE REGRESIÓN trControl = control <- trainControl(method='repeatedcv', number=10, search = 'grid')) Para la variable capturas se realizó una agrupación por la variable Fecha que engloba Mes y Año, se calculó el número de capturas por trampas ocurridas en cada localización y mes. Posteriormente se halló la media por cada localización, y se realizó la ponderación para quedarse con la variable “capt” que lo que expresa es si el número de capturas fue superior o inferior a la media y en qué medida. Aquellas localizaciones con menos de tres registros a lo largo de todo el estudio fueron descartadas. Por último, se unieron con el conjunto de los registros de variables ambientales, en función de los que compartían fecha y localización. ev_pres_abs <- ev %> % mutate(Fecha=year(eventDate)+(month(eventDate)-1)/12, IDL=localityID) %> % group_by(IDL, Fecha) %> % summarise(tot_capt = sum(sampleSizeValue,na.rm=TRUE)) %> % ungroup() %> % group_by(IDL) %> % mutate(med = mean(tot_capt, na.rm=TRUE), capt = tot_capt/med, cont= n()) %> % filter(cont > 3) %> % select(-tot_capt, -med, -cont) %> % ungroup() %> % inner_join(pres_abs, by=c("IDL","Fecha")) Los datos anteriores tienen un tamaño considerable, aunque en el caso de las variables ambientales el porcentaje de registros incompletos era superior a la mitad. Para solucionar esto, al aplicar el preprocesado de los datos mediante recipe se incluyó la imputación de las variables predictoras en función del resto de predictores numéricos. Aparte de eso, se eliminaron aquellas variables con una correlación superior a 0.8 en valor absoluto y se realizó un centrado y escalado de los predictores numéricos. lagunas_imputado <- ev_pres_abs %> % filter(Habitat == "Lagunas temporales") %> % select(-Habitat, -Salinidad, -Fecha, -IDL) %> % recipe(capt ~ .) %> % step_corr(all_numeric(), threshold = 0.8) %> % step_normalize(all_numeric_predictors()) %> % step_impute_knn(all_predictors(), impute_with = imp_vars(all_numeric_predictors(), -capt)) %> % prep() %> % bake(new_data = NULL) Se realizó este mismo procedimiento para las marismas. De igual forma, para ambos hábitats se aplicó el preprocesado a los datos pero sin realizar la imputación, de modo que pudiera compararse el rendimiento de los modelos creados en base a cada uno de los conjuntos. 60 CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 5.4. MODELOS DE REGRESIÓN La división en entrenamiento y test se llevó a cabo igual que en los casos anteriores, utilizando la división 80/20 también en estos casos. En este caso se trabajó con el paquete workflows para generar los modelos de bosques aleatorios y de k-vecinos más cercanos, ambos con su correspondiente ajuste de hiperparámetros. Se comenzó creando los modelos base, indicando los parámetros que serán ajustados: rf_lag_capt_tune_model <- rand_forest(mtry = tune(), trees = tune(), min_n = tune()) %> % set_engine("randomForest") %> % set_mode("regression") knn_tune_model_reg <- nearest_neighbor(neighbors = tune(), weight_func = tune(), dist_power = tune()) %> % set_engine("kknn") %> % set_mode("regression") Este modelo se añade a un workflow en blanco, así como la fórmula indicando cuál es la variable de interés. Para el ajuste de hiperparámetros se utilizan valores aleatorios gracias a grid_random, y se aplica validación cruzada con tune_grid usando los pliegues proporcionados por vfold_cv. Por último con select_best se indica la medida que se quiere utilizar para decidir cuál es la mejor combinación de parámetros, se aplica finalize_workflow y con last_fit estaría listo para obtener las métricas y las predicciones. El código ha sido explicado de forma reducida, pero puede encontrarse por completo el código de creación de todos los modelos en el Apéndice A. CAPÍTULO 5. INTRODUCCIÓN A LOS MODELOS PREDICTIVOS 61 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.10: Matriz de confusión regresión logística, sexo vs largo cuerpo, proporción Female Male Female 50.95 49.05 Male 44.13 55.87 Tabla 6.11: Matriz de confusión regresión logística, sexo vs variables tamaño Female Male Female 64 63 Male 48 104 y Hembras en esta especie no radican en el tamaño de los mismos. En caso de que sea necesario poder clasificar los individuos sin tener que examinar especificamente el sexo de cada uno, este estudio nos indica que debería barajarse qué otras variables podrían servir para este cometido, ya que las disponibles en este conjunto no son suficientes. Se pueden reflejar las medidas de rendimiento de cada uno de los modelos para compararlos. Concretamente en la Tabla 6.13 se muestra el modelo con todas las variables de tamaño de cada tipo. En las tres medidas el modelos de bosques aleatorios tiene mejores resultados. Aparte, la sensibilidad del modelo de Regresión Logística es muy mala, por lo que podría decirse que este modelo es el peor aunque tenga mejor especificidad que el de K vecinos más cercanos. Ninguno de los tres modelos obtiene medidas como para poder sacar conclusiones fiables de ellos. 6.2. Predicción de la Madurez El proceso con esta variable es el mismo que se siguió con la variable Sexo, primero se divide de aquellos individuos que se conoce su variable Madurez y que no tenían ningún valor perdido en las variables de tamaño en dos muestras: entrenamiento y test. Lo mismo se realiza sobre el conjunto de aquellos para los que se conoce la variable Madurez, pues de todos ellos se dispone de la variable Largo del cuerpo, y por tanto se puede intentar realizar la predicción solo con esa variable y teniendo un tamaño muestral mucho mayor. En este caso se dispone de información sobre 8968 individuos inmaduros, 6503 maduros y 11 sin determinar. Por tanto, al restringir aquellos de los que sí se dispone información y dividir en dos muestras se obtendrían: 6277 individuos inmaduros y 4552 individuos maduros en la muestra de entrenamiento, 2691 individuos inmaduros y 1951 maduros en la muestra test. Por otro lado, una vez se obtienen solo los que tienen registradas sus cuatro variables de tamaño, se reduce a: 447 individuos inmaduros y 542 maduros. Se han invertido los tamaños entre ambas categorías. Por tanto, para la muestra de entrenamiento la estructura Tabla 6.12: Matriz de confusión regresión logística, sexo vs variables tamaño, proporción Female Male Female 50.39 49.61 Male 31.58 68.42 68 CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.13: Comparativa modelos determinación Sexo Modelo Exactitud Especificidad Sensibilidad K vecinos más cercanos 0.645 0.678 0.614 Bosques aleatorios 0.685 0.704 0.677 Regresión logística 0.602 0.684 0.504 Tabla 6.14: Matriz de confusión KNN, Madurez vs variables de tamaño Immature Mature Immature 120 15 Mature 24 139 es 312 inmaduros y 379 maduros, y para la muestra de test es de 135 inmaduros y 163 maduros. 6.2.1. K vecinos más cercanos Se debe determinar el k óptimo para cada uno de los dos modelos, para ello se utiliza la muestra de entrenamiento para comprobar los valores impares que van de 5 a 25 ambos inclusive. Se obtiene el k óptimo para el modelo con las cuatro variables. ## [1] "El mejor k para el conjunto de entrenamiento es 17" Una vez se tiene ese parámetro se pueden realizar las predicciones para la muestra test y comparar con los valores reales para ver su precisión. La precisión de este modelo está alrededor del 90% en la muestra test, así que cabe esperar resultados similares si se utilizase para clasificar nuevos individuos en el futuro. Es un modelo mucho más fiable que el que se creó para clasificar Machos y Hembras, pero antes de sacar conclusiones se comprobarán también los demás métodos antes utilizados. Para el modelo utilizando solo la variable Largo del cuerpo también se obtiene el número óptimo de vecinos a utilizar. ## [1] "El mejor k para el conjunto de entrenamiento es 21" Utilizando este parámetro se predice la muestra test, y los niveles de acierto resultantes son los que se indican en las Tablas 6.16 y 6.17. Incluso usando solo una variable para predecir el estado de Madurez de los cangrejos, los resultados tienen una fiabilidad del 90% aproximadamente. Este modelo y el anterior pueden ser considerados fiables y podrían utilizarse en un futuro en caso de tener nuevos datos que se quieran registrar, o para clasificar aquellos valores perdidos que había en el conjunto disponible. Tabla 6.15: Matriz de confusión KNN, Madurez vs variables de tamaño. Proporción Immature Mature Immature 88.89 11.11 Mature 14.72 85.28 CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 69 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.16: Matriz de confusión KNN, Madurez vs largo cuerpo Immature Mature Immature 2467 224 Mature 192 1759 Tabla 6.17: Matriz de confusión KNN, Madurez vs largo cuerpo. Proporción Immature Mature Immature 91.68 8.32 Mature 9.84 90.16 6.2.2. Bosques aleatorios y árboles de clasificación En este caso, al igual que para clasificar Hembras y Machos, se presenta un ejemplo de un árbol de clasificación utilizando únicamente la variable Largo del cuerpo, y además se crea un modelo mediante bosques aleatorios utilizando las cuatro variables. largo_cuerpo 1 ≤68 >68 Node 2 (n = 6022) Mature Immature 0 0.2 0.4 0.6 0.8 1Node 3 (n = 4807) Mature Immature 0 0.2 0.4 0.6 0.8 1 Hay un único nodo en este caso, esto se debe a que la distribución era mucho más clara para esta variable: los individuos inmaduros presentan generalmente un menor tamaño y los maduros un tamaño mayor. Los nodos finales presentan una distribución mucho menos balanceada que en el caso de la variable sexo, por lo que al menos para la muestra de entrenamiento se obtiene una clasificación aparentemente buena. El punto clave estaba en fijar un punto crítico que permita hacer una buena división de los grupos. Al aplicar esto al conjunto entrenamiento y al conjunto test se obtienen los resultados mostrados en las Tablas 6.18 y 6.19. 70 CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.18: Matriz de confusión árbol clasificación, Madurez vs largo, entrenamiento Immature Mature Immature 90.25 9.75 Mature 7.84 92.16 Tabla 6.19: Matriz de confusión árbol clasificación, Madurez vs largo, test Immature Mature Immature 90.52 9.48 Mature 8.20 91.80 En ambos casos es muy similar la tasa de aciertos, lo cuál confirma que no se trata de un modelo sobreajustado, sino que su fiabilidad es verídica. Tras haber visto esto, para el conjunto de todas las variables de tamaño se utiliza el método de bosques aleatorios. Los resultados obtenidos son los de las Tablas 6.20 y 6.21. El nivel de fiabilidad que tiene este modelo es similar al del árbol antes creado, ambos están en torno al 85-90% de acierto. En caso de tener todas las variables es conveniente utilizar este modelo ya que tiene mayor información, pero si solo se dispone del largo del cuerpo podría utilizarse como norma de clasificación la establecida en el árbol. 6.2.3. Regresión lógística Finalmente, se plantea diseñar un modelo de regresión logística, pudiendo comparar la eficacia que tenga el modelo que utilice solo la variable Largo del cuerpo respecto al que las usa todas, así como comparar la eficacia de estos con sus equivalentes en las técnicas anteriores. La predicción utilizando solo la variable Largo del cuerpo da mejores resultados para los individuos inmaduros, aunque en ambas categorías está alrededor del 90%, como era de esperar viendo los modelos anteriores. Para el modelo de regresión logística con todas las variables tamaño como predictoras, los resultados obtenidos se muestran en las Tablas 6.24 y 6.25. Este modelo tiene un rendimiento ligeramente peor que el de una sola variable, pero puede deberse a la selección de la muestra. Lo que no es coincidencia es que los modelos que predicen la etapa de Madurez de los cangrejos sean mucho más fiables que los que predicen el Sexo. Al menos los basados en las variables de tamaño, que son los que se han podido estudiar con la información proporcionada. Por tanto, se concluye que en caso de quererse predecir en un futuro la etapa de Madurez de los cangrejos, bastaría con medir las variables de tamaño, o incluso solo el Largo Tabla 6.20: Bosque clasificación, matriz de confusión Immature Mature Immature 119 16 Mature 24 139 CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 71 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.21: Bosque clasificación, matriz de confusión proporción Immature Mature Immature 88.15 11.85 Mature 14.72 85.28 Tabla 6.22: Matriz de confusión regresión logística, Madurez vs largo cuerpo Immature Mature Immature 2492 199 Mature 223 1728 Tabla 6.23: Matriz de confusión regresión logística, Madurez vs largo cuerpo, proporción Immature Mature Immature 92.60 7.40 Mature 11.43 88.57 Tabla 6.24: Matriz de confusión regresión logística, Madurez vs variables tamaño Immature Mature Immature 122 13 Mature 24 139 Tabla 6.25: Matriz de confusión regresión logística, Madurez vs variables tamaño, proporción Immature Mature Immature 90.37 9.63 Mature 14.72 85.28 72 CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 6.2. PREDICCIÓN DE LA MADUREZ Tabla 6.26: Comparativa modelos determinación Madurez Modelo Exactitud Especificidad Sensibilidad K vecinos más cercanos 0.909 0.9 0.916 Árbol de clasificación 0.911 0.918 0.905 Regresión logística 0.909 0.886 0.926 del cuerpo, para poder conseguir una estimación fiable. En cambio, si el objetivo del estudio fuera poder predecir si el cangrejo en cuestión es Macho o Hembra, se necesitaría información adicional que estuviera relacionada con otras características de la especie. Esta cuestión podría trasladarse a los ecólogos en caso de que se diera otra colaboración en el futuro, y tal vez replantear las mediciones tomadas o añadir alguna más para facilitar la predicción de la variable Sexo. También sería interesante saber si el hecho de que la mayoría de valores perdidos para esta variable se encuentren en individuos Inmaduros es por la dificultad de identificar a estos cuando aún son demasiado pequeños. Aunque el modelo matemático reduzca estas cuestiones a números y porcentajes, no hay que perder de vista que se trata de una situación real y que involucra distintas disciplinas, de modo que intentar obtener resultados desde un punto de vista puramente matemático puede provocar que el enfoque no sea el más eficaz para la situación. Por último, para resumir esta sección se refleja el mejor modelo de cada tipo en una tabla comparativa (Tabla 6.26). En este caso, el mejor ajuste se da utilizando únicamente la variable “Largo cuerpo”. El árbol de clasificación tiene la mejor Exactitud y Especificidad, pero la peor Sensibilidad. En conjunto parece que los tres modelos tienen resultados similares, por lo que cualquiera de ellos sería válido. CAPÍTULO 6. MODELOS DE CLASIFICACIÓN 73 Capítulo 7 Modelos de regresión Tratar con las variables numéricas supone una gran diferencia con respecto a las variables categóricas. Los métodos antes utilizados podrían servir, pudiendo cambiar la regresión logística por regresión lineal. Sin embargo, la principal diferencia se encuentra en la información disponible y la forma de tratarla y combinar los diferentes conjuntos. En el caso anterior se intentaba predecir información sobre un cangrejo concreto dados unos valores de ese mismo cangrejo (longitud, peso, ancho). El primer apartado de este capítulo trata de eso mismo: intentar predecir una de las variables de tamaño de un cangrejo en función del resto. Pero no siempre ocurre esto, especialmente cuando se trabaja con el conjunto de datos ambientales. Es importante recordar que este conjunto de datos no ha sido publicado al mismo tiempo que los otros tres que se manejan en este estudio, y por tanto es un fichero con una estructura diferente y un trabajo de pre-procesado bastante menor por parte de las personas que recopilaron los datos. Gran parte del esfuerzo dedicado a este capítulo radicó en intentar pulir el conjunto de datos y aprovechar al máximo la información disponible, aun cuando se encontraron las dificultades que a continuación se explican. En primer lugar, los datos publicados trataban con fechas en días concretos, ya fuera el día de la captura de un cangrejo concreto o el día en que se revisaron las trampas; todos estos valores estaban asignados a una fecha. En cambio, en el fichero con los datos ambientales los datos se registran por meses: se indica una localización, un año y un mes, y se registran los valores de las variables, pero no queda registrado en qué día en concreto se realizó la medición. La forma de solucionar esto fue agrupar las otras variables por meses también, de modo que ambos conjuntos pudieran combinarse en los casos que coincida la fecha y la localización. En segundo lugar, aunque en cada mes y localización se tomase un único dato de cada variable (en caso de haberse tomado alguno), examinando el fichero puede verse que hay un gran número de registros que contienen la misma información. La única diferencia entre estos registros, tal como se explicó en el apartado de lectura de datos, era el código de nasa que se asignaba, pero eso no influye en las variables ambientales y por tanto es una información que fue descartada. Había también casos en los que, para un mismo mes y localización, se registraban datos prácticamente iguales, pero que diferían en una única variable. Estos datos fueron considerados erróneos y eliminados por completo, ya que daban lugar a incongruencias, pudiendo llegar a tomar valores simplemente imposibles a nivel científico como era el caso del pH. En general, esta limpieza convirtió el conjunto original de 6710 registros en uno de 784, pero ya sin repeticiones ni ruido. 75 7.1. TAMAÑO DE LOS CANGREJOS En tercer lugar, de estos registros solo 197 eran completos. Todos los demás tenían algún valor perdido, en la mayoría de casos una cantidad considerable de ellos. En el caso de la variable Nitrato, si se eliminase se conseguiría llegar a 221 registros completos, pero eso solo ocurre ya que es la que mayor porcentaje de valores perdidos tenía, y no funcionaría intentar seguir aplicando una estrategia similar. El problema que esto supone es que a la hora de hacer modelos se deben separar los datos de cada tipo de hábitat, ya que tienen comportamientos muy distintos tanto en valores de estas variables como en las capturas en trampas. Haciendo esto los conjuntos que quedan se componen de 111 registros completos para las Lagunas y 86 para las Marismas. Si a eso se le añade que debe hacerse una división en conjunto de entrenamiento y conjunto de prueba (o test) para poder evaluar el rendimiento de los modelos, esto hace que la cantidad de datos con los que se podía trabajar fuera insuficiente. Una forma de solucionarlo podría ser aplicar validación cruzada a los modelos, pero los valores perdidos siguen suponiendo un enorme desperdicio de información ya que en modelos como árboles aleatorios aquellas observaciones que no están completas son descartadas directamente. Este último problema no tiene una solución como tal, pero la forma de intentar remediarlo que se adoptó en este estudio fue imputar los datos perdidos mediante el método de los k-vecinos más cercanos como se explicó en el capítulo introductorio de este bloque. Gracias al paquete recipes del ecosistema tidymodels de R, esta tarea se llevó a cabo de forma automatizada, y los conjuntos finales tienen 143 registros en el caso de las Lagunas y 140 registros en el caso de las Marismas, lo cual sigue siendo bastante menos de lo que sería ideal, pero al menos es un tamaño considerable. A pesar de todo, a la hora de establecer una relación entre las variables ambientales y las capturas realizadas no está claro que vaya a darse una buena modelización limitándose a relacionar aquellos datos coetáneos. Tal vez el número de cangrejos se vea más influenciado por un acumulado de estas variables en los meses (o semanas) anteriores, y con los datos disponibles no hay forma de llevar a cabo ese estudio. En la mayoría de localizaciones los datos no se corresponden con meses consecutivos, sino que se trata de registros dispersos a lo largo del tiempo en que se llevó a cabo la recogida de datos, habiendo a lo sumo dos meses consecutivos en aquellas localizaciones en las que se registraron mayor número de registros. Por tanto, para poder cuadrar los datos de la forma que se consideró más adecuada, se combinaron los datos que coincidían en localización y fecha, y se intentó obtener modelos apropiados en base a esto. 7.1. Tamaño de los cangrejos De las cuatro variables de tamaño sobre los cangrejos, una de ellas (largo del cuerpo) se midió para todos los individuos; en cambio, las otras tres apenas se midieron en el 10%. Concretamente, el largo del cefalotórax se midió en el 12.5% de los cangrejos que se estudiaron, y el ancho de cefalotórax y el peso en el 7%. Sería de interés por tanto ver si es posible estimar alguna de estas variables respecto al largo del cuerpo. Para cada una de las tres variables se trabaja sobre el subconjunto de cangrejos en los que esta variable fue medida, aplicando un modelo de regresión lineal solo con una sola variable predictora: el largo del cuerpo. El error cometido se calcula mediante Validación Cruzada (fijando una semilla) de modo que se mida sobre la muestra de entrenamiento, y la muestra test se utilice para visualizar las predicciones y compararlas con los valores 76 CAPÍTULO 7. MODELOS DE REGRESIÓN 7.1. TAMAÑO DE LOS CANGREJOS reales. Ya que cada una de las tres variables tiene rangos distintos, comparar el error en términos absolutos podría ser una forma engañosa de comparar los modelos. El método elegido para evaluar el error cometido por los modelos es calcular el RMSE por Validación Cruzada, y calcular el porcentaje que este error supone con respecto al valor medio de la variable en el conjunto de entrenamiento. Antes de empezar a hacer ningún tipo de modelización, es conveniente analizar las correlaciones que tienen las variables. Utilizando para cada par de variables todos los registros completos de ambas, la matriz de correlaciones obtenida se ve en la Figura 7.1. −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1 largo_cuerpo largo_cefalotorax ancho_cefalotorax peso largo_cuerpo largo_cefalotorax ancho_cefalotorax peso Figura 7.1: Correlaciones variables tamaño A primera vista llama la atención que el largo del cefalotórax y el ancho del cefalotórax sean el par que menos correlación tienen entre sí. Concretamente, su correlación es -0.0698 mientras que la siguiente menor (en valor absoluto) es 0.515 entre largo del cefalotórax y peso. El motivo de que esto ocurra no tiene explicación biológica, pero mirando a fondo los datos puede observarse un fenómeno que se refleja en la Figura 7.2. 0 20 40 60 80 2005 2010 2015 fecha largo_cefalotorax Largo del cefalotórax vs fecha 0 20 40 60 2008 2010 2012 2014 2016 fecha ancho_cefalotorax Ancho del cefalotórax vs fecha Figura 7.2: Tamaños cefalotórax según fecha Los datos parecen dividirse en dos periodos. Por poner un punto de separación, podría distinguirse el grupo de cangrejos cuyas medidas se tomaron antes de 2010 y el grupo de CAPÍTULO 7. MODELOS DE REGRESIÓN 77 7.1. TAMAÑO DE LOS CANGREJOS el eje y la predicción; eso no significa que la variable objetivo tenga esa distribución. Aun así, es razonable plantearse si al hacer alguna transformación logarítimica en el modelo anterior los resultados obtenidos serían mejores; por tanto, se procede a estudiar dicho modelo. ## [1] "El error cometido con la transformación logarítmica es del 16.19 %" El error cometido en la predicción del logaritmo del peso (respecto a la media del logaritmo del peso en la muestra test) es mucho menor que el que se cometía en la predicción del peso normal. Eso es buena señal, ya que los resultados concuerdan con las gráficas. 0 25 50 75 100 0 20 40 60 peso pred Predicción del peso (destransformado) frente a valor real 0 25 50 75 100 0 20 40 60 peso pred Filtrando los casos con predicciones < 100 Predicción del peso (destransformado) frente a valor real Figura 7.10: Predicciones regresión con transformación logarítmica peso Para la visualización se enfrenta la predicción del peso (la exponencial de la estimación del logaritmo que se obtuvo mediante la regresión lineal) respecto al valor real del peso. Al haber tenido que aplicar la exponencial, hay algunos casos cuya predicción se dispara, y hacen que la gráfica quede distorsionada. Para poder visualizar mejor la parte de la gráfica que concentra la mayor parte de los casos, en la segunda gráfica se observan solo aquellos que presentaron una predicción inferior a 100 gr. En ambas puede verse que, aunque algunos casos se alejen bastante del valor real, en la mayoría de los individuos la predicción mejora de forma significativa. 7.1.2. Bosques aleatorios Otra alternativa al modelo inicial es intentar utilizar un modelo de regresión distinto, como el método de bosques aleatorios, para poder así utilizar las tres variables predictoras en vez de solo una. Para este modelo se realiza el ajuste sobre el parámetro “mtry”, que indica cuántas variables deben estudiarse en cada nodo para decidir la mejor partición, obteniendo en este caso como resultado que el valor óptimo es 1. Por tanto, en cada nodo se elige aleatoriamente una de las tres variables predictoras, y se divide de la mejor forma posible usando esta variable. Para el número de árboles generados en el bosque se deja el valor por defecto, ya que con el tamaño del conjunto de datos disponible tampoco es necesario generar un enorme número de árboles, y la mejora que podría conseguirse sería a conllevaría un aumento del coste computacional excesivo. 84 CAPÍTULO 7. MODELOS DE REGRESIÓN 7.1. TAMAÑO DE LOS CANGREJOS Con estos parámetros, se calcula el porcentaje de error cometido sobre la muestra entrenamiento, al igual que se hizo con los modelos anteriores, y se muestra la gráfica de las predicciones de la muestra test frente a los valores reales. ## [1] "El error cometido es del 14.34 %" 0 20 40 60 0 30 60 90 120 peso pred Figura 7.11: Predicciones bosque aleatorio: peso según variables tamaño Puede verse que el error cometido con el modelo de bosques aleatorios es mucho me- nor que con el modelo de regresión, en parte porque han podido aprovecharse las tres variables disponibles, y también porque este método tiene un proceso más complejo de entrenamiento que permite hacer un ajuste más preciso. También se puede ver gráficamente cómo las predicciones se acercan más a la línea diagonal que representa el ajuste perfecto. De este apartado pueden sacarse varias conclusiones: Alrededor de la pausa que hubo en las mediciones de largo del cefalotórax y en ancho del cefalotórax de los cangrejos en 2010 los comportamientos de estas variables se vieron modificados de manera considerable. En cambio esto no ocurrió para la variable peso, que también tuvo un parón en las mediciones tomadas, ni para la variable largo del cuerpo, que fue tomada de forma continuada durante los años del estudio. El motivo de este suceso es desconocido. Es conveniente, para realizar estimaciones de las variables relacionadas con el cefalotórax, tener en cuenta en cuál de los dos periodos se encuentra el individuo de interés. Para predicciones futuras deberían utilizarse únicamente los valores posteriores a 2010 como conjunto de entrenamiento, pues los datos anteriores producen una perturbación que empeora la calidad de las predicciones. Al estimar la variable peso mediante un modelo de regresión con la variable largo del cuerpo se comete un error muy superior a lo que sería adecuado. En caso de CAPÍTULO 7. MODELOS DE REGRESIÓN 85 7.2. NÚMERO DE CAPTURAS EN TRAMPAS Tabla 7.6: Comparativa modelos regresión lineal según periodo Modelo R^2 ajustado % error RMSE Largo Cefalotórax, antes de 2010 0.86514 19.22 2.768 Largo cefalotórax, después de 2010 0.94821 8.97 3.051 Ancho Cefalotórax, antes de 2010 0.78419 12.31 5.11 Ancho cefalotórax, después de 2010 0.93693 10.58 1.586 Tabla 7.7: Comparativa modelos predicción peso Modelo % error RMSE Regresión logarítmica, Peso 16.19 1.444 Bosque aleatorio, Peso 14.34 2.403 querer obtener predicciones de esta variable, es mejor aplicar una transformación logarítmica antes de aplicar el modelo de regresión lineal. Por otro lado, también se obtienen buenas predicciones utilizando el método de bosques aleatorios, que genera un nivel de error similar. En el caso de la regresión lineal, pueden darse valores puntuales más alejados de lo deseado, pero solo requiere de la variable largo del cuerpo como predictora. En cambio, los bosques aleatorios proporcionan una estimación menos propensa a generar predicciones alejadas en extremo, pero requiere de las tres variables de tamaño para obtener la predicción. A modo de resumen se muestran los resultados de los mejores modelos para las variables Largo de Cefalotórax y Ancho del Cefalotorax (ver Tabla 7.6). Y lo mismo para los modelos del peso en Tabla 7.7. Ya que el porcentaje de error depende de la media de la variable para el conjunto en que fue estudiada (entrenamiento), es una medición del error disinta al RMSE, por tanto se muestran ambos para cada modelo realizado. 7.2. Número de capturas en trampas Para poder utilizar las variables ambientales, como se explicaba al inicio del capítulo, lo que debe hacerse es agrupar los registros diarios en registros mensuales, de forma que puedan coincidir con los datos de las variables ambientales que fueron tomados. También se filtraron aquellas localizaciones en las que a lo largo de todo el estudio tuvieron 3 o menos registros, ya que una localización en la que se haya puesto una trampa o dos no aporta información suficiente como para ser tenida en cuenta. La variable objetivo son las capturas ponderadas, obtenidas dividiendo el número real de capturas entre la media de la localización correspondiente. Para estudiar esta variable respecto a las variables climáticas, se obtienen los valores ponderados de forma mensual por cada localización. Una vez realizadas estas transformaciones, se combinan ambos conjuntos y se realiza la imputación. Al mismo tiempo, se normalizan las variables numéricas predictoras y se eliminan aquellas que tienen demasiada correlación. En el caso de la Salinidad y la Conductividad, dada la enorme cantidad de datos perdidos de la primera se elimina esta de forma manual para evitar perder información innecesariamente. Tras obtener el conjunto de datos ambientales sin valores perdidos, puede combinarse con el fichero de capturas en trampas y proceder a la modelización. 86 CAPÍTULO 7. MODELOS DE REGRESIÓN 7.2. NÚMERO DE CAPTURAS EN TRAMPAS 7.2.1. Bosques aleatorios El número de variables predictoras en este caso es muy superior al del apartado anterior. Incluso eliminando aquellas con correlaciones demasiado altas, siguen quedando 13 variables, por lo cual el método de bosques aleatorios es bastante conveniente, ya que solo utiliza un subconjunto de estas en cada división. Los parámetros que pueden optimizarse son: mtry (número de variables elegidas para cada nodo), trees (número de árboles en el bosque) y min_n (tamaño mínimo para que un nodo vuelva a dividirse). El tuning de estos parámetros se realiza mediante tidymodels. Se evalúa tanto para los datos de las lagunas temporales como para las marismas el porcentaje de error cometido sobre la variable de capturas ponderadas en la muestra test. ## [1] "El valor óptimo de mtry en las lagunas es 1" ## [1] "El error cometido en las lagunas es del 25.84 %" ## [1] "El valor óptimo de mtry en las marismas es 1" ## [1] "El error cometido en las marismas es del 27.55 %" En ambos casos el error está alrededor del 20 %. Es bastante grande, pero hay que considerar que el valor medio de esta variable es 1, de modo que cualquier pequeña interferencia supone un porcentaje de error considerable. 0.5 1.0 1.5 0.5 1.0 1.5 Numero capturas real Predicción numero capturas Lagunas 0.5 1.0 1.5 0.5 1.0 1.5 Numero capturas real Predicción numero capturas Marismas Figura 7.12: Predicciones capturas bosques aleatorios según hábitat En las gráficas se puede ver la comparación entre los valores reales de la variable y los valores de sus predicciones. Para poder hacer una comparación adecuada, ambas se han representado con una misma escala. En ninguna de las dos gráficas hay demasiados individuos sobre la diagonal, y no se ve que los comportamientos sean muy diferentes entre ambos modelos. El caso de las lagunas tiene casos más dispersos pero el resto están más concentrados cerca de la diagonal, lo cual compensa el error cometido y por ello ambos modelos tenían un porcentaje de error similar. Estos dos modelos fueron realizados sobre datos de las variables predictoras que en algunos casos fueron imputados, esto favorece una muestra de mayor tamaño pero tal vez CAPÍTULO 7. MODELOS DE REGRESIÓN 87 7.2. NÚMERO DE CAPTURAS EN TRAMPAS provoque una estimación incorrecta. En caso de haber tenido suficientes datos al descartar aquellos con valores perdidos, hubiera sido preferible limitarse a usar estos. Para confirmar que la decisión de imputar los valores perdidos de las variables predictoras fue correcta, se generan también los correspondientes modelos utilizando únicamente los registros completos. El error cometido en ese caso es: ## [1] "El valor óptimo de mtry en las lagunas es 1" ## [1] "El error cometido en las lagunas es del 12.81 %" ## [1] "El valor óptimo de mtry en las marismas es 1" ## [1] "El error cometido en las marismas es del 11.41 %" El error cometido parece ser menor, aunque cabe destacar que en el proceso de creación y evaluación del modelo fueron numerosas las alertas (warnings) que saltaron, advirtiendo de que se requerían muestras de mayor tamaño que las disponibles. Esto significa que el modelo generado no cumple los mínimos con los que fueron ideados los comandos de R utilizados, y por tanto no tienen garantías de ser correctos. 0.8 0.9 1.0 1.1 1.2 1.3 0.8 0.9 1.0 1.1 1.2 1.3 Numero capturas real Predicción numero capturas Lagunas 0.8 0.9 1.0 1.1 1.2 1.3 0.8 0.9 1.0 1.1 1.2 1.3 Numero capturas real Predicción numero capturas Marismas Figura 7.13: Predicciones capturas bosques aleatorios según hábitat Visualmente lo primero que se aprecia es la cantidad de puntos en cada una de las gráficas, mucho menor que en el anterior caso. Además, en el caso de las lagunas casi todas las predicciones fueron iguales, incluso cuando el rango de valores reales era amplio para la variable con la que se está tratando. En el caso de las marismas esto no es tan llamativo, pero también se nota que los datos están más dispersos que en el modelo creado en base a los datos imputados. Por tanto, el uso de datos imputados para las variables predictoras que tenían valores perdidos está justificado en este caso, y no haber recurrido a este método podría haber llevado a error con modelos que engañosamente parecen tener un error pequeño pero que no proporcionan estimaciones de calidad. 88 CAPÍTULO 7. MODELOS DE REGRESIÓN 7.2. NÚMERO DE CAPTURAS EN TRAMPAS 7.2.2. K vecinos más cercanos Por último, puede implementarse el método de los K vecinos más cercanos. De nuevo se realizarán los modelos de forma separada en cada uno de los hábitats y partiendo del conjunto de datos con los valores perdidos de las variables predictoras imputados. No se crea en este caso el modelo con los datos sin imputar, ya que se vio que no es correcto utilizar ese conjunto. Se realiza tuning para los parámetros, de modo que se escojan los mejores para la muestra de entrenamiento generada, y el porcentaje de error se estudia sobre la muestra test. Para el conjunto de datos de las lagunas los resultados obtenidos son: ## [1] "El mejor k para el conjunto de entrenamiento de las lagunas es 9" ## [1] "El error cometido en las lagunas es del 25.51 %" Y lo mismo para el conjunto de las marismas: ## [1] "El mejor k para el conjunto de entrenamiento de las marismas es 14" ## [1] "El error cometido en las marismas es del 26.48 %" El error cometido en porcentaje es cercano al que se cometía con los modelos mediante bosques aleatorios, no hay una diferencia significativa como para afirmar que uno de los métodos sea mejor que el otro. 0.5 1.0 1.5 0.5 1.0 1.5 Capturas real Predicción capturas Lagunas 0.5 1.0 1.5 0.5 1.0 1.5 Capturas real Predicción capturas Marismas Figura 7.14: Predicciones capturas KNN según hábitats Las predicciones parecen tener una distribución similar a las obtenidas en los modelos de bosques aleatorios en ambos casos. Para las lagunas están más agrupadas en la zona central de la diagonal pero con algún caso dispreso. Para las marismas están más repartidas en general, aunque el error medio cometido es casi igual en ambos conjuntos. De todos modos, como ya se comentó al inicio del capítulo, han sido muchos los problemas afrontados para conseguir obtener modelos coherentes dados los conjuntos de datos disponibles. Probablemente si se tuvieran datos de forma más continuada podrían haberse CAPÍTULO 7. MODELOS DE REGRESIÓN 89 7.2. NÚMERO DE CAPTURAS EN TRAMPAS Tabla 7.8: Comparativa modelos sobre capturas Modelo % error RMSE Bosques aleatorios, lagunas 25.843 0.267 Bosques aleatorios, marismas 27.553 0.278 K vecinos más cercanos, lagunas 25.506 0.263 K vecinos más cercanos, marismas 26.482 0.268 implementado otros métodos que permitieran relacionar la cantidad de capturas de un mes en concreto no solo con las variables registradas ese mes sino con los meses anteriores. No se tienen garantías de que ese modelo fuera a cometer un error menor que los aquí presentados, pero de cara al significado biológico de los datos tendría una mayor sentido. Al igual que en la sección anterior, los mejores modelos se representan en la Tabla 7.8. Hay dos modelos por cada hábitat, solo se han mostrado aquellos que utilizaron los datos imputados climáticos. Los resultados son muy similares para cada par de modelos, de modo que podrían utilizarse indistintamente para obtener estimaciones según el hábitat de interés. 90 CAPÍTULO 7. MODELOS DE REGRESIÓN Parte IV Conclusiones 91 Capítulo 8 Conclusiones En este último capítulo del documento se pretende reflejar cómo fueron abordados los objetivos iniciales presentados en la Introducción, los pasos que se dieron a lo largo del estudio y las conclusiones que fueron extraídas. Antes de entrar en detalle sobre los modelos obtenidos, se hará un repaso de aquellos puntos de interés aportados en este trabajo como medios para alcanzar dichos objetivos. También se mencionarán las vías de trabajo que no pudieron llevarse a cabo, y aquellas que en un trabajo posterior podría plantearse emprender. 8.1. Aportaciones iniciales Para el estudio inicial, se disponía de cuatro conjuntos de datos. Tres de ellos con un formato más depurado que fueron publicados simultáneamente, y otro menos preparado que requería de un tratamiento previo para poder ser estudiado. La información disponible trataba sobre capturas individuales de cangrejos, así como de capturas en trampas e información sobre el clima y las variables ambientales en las disitntas localizaciones de Doñana que fueron estudiadas. Estos conjuntos fueron examinados en detalle, eliminando aquella información redundante o errónea, y se realizaron múltiples transformaciones para poder trabajar mejor con ellos. Esta limpieza y ordenación de los datos permitió tener conjuntos más reducidos pero ya listos para ser tratados, en caso de querer realizar un estudio con un enfoque distinto. Es importante destacar la creación de una nueva variable: capturas ponderadas. Esta variable fue creada para permitir tratar con la variable de capturas en trampas de una forma distinta, prestando mayor atención a si hubo más o menos capturas de lo habitual en cada sitio y no solo al número de capturas, ya que esta variable no es en sí un fin sino un medio para estimar la población total de cangrejos y su evolución, cosa que difícilmente se puede estimar a partir de los datos disponibles. Al haber enfrentado un problema real, esto conllevó tratar con problemas en la metodología y continuidad de la recopilación de datos, además de la considerable cantidad de datos perdidos. Es debido a esto que el estudio a través de Series Temporales, entre otros, tuvo que ser descartado. Fue necesario revisar en detalle, importar y tratar cuidadosamente los conjuntos, relacionándolos entre sí cuando fuera posible, o reestructurándolos cuando fuera necesario. Se estudiaron a fondo las variables, para discernir cuáles aportaban información y cuales no eran de provecho para el estudio llevado a cabo. Para 93 A.1. IMPORTACIÓN mof_bruto <- read_delim("data/Don_ biom_red-swamp-crayfish _mof_20230111.csv", delim = ";",escape_double = FALSE, trim_ws = TRUE) mof <- mof_bruto %> % select(occurrenceID, measurementValue, measurementType) %> % separate(occurrenceID, sep="_", into=c("IDL","fecha","Cangrejo")) %> % mutate(fecha = as.Date(fecha, format=" %Y- %m- %d")) %> % pivot_wider(names_from = measurementType, values_from = measurementValue) %> % rename(largo_cuerpo = `total body length`, largo_cefalotorax = `cephalothorax length`, ancho_cefalotorax = `cephalothorax width`, peso = weight) %> % filter(largo_cuerpo < 900) %> % mutate(IDL = as.numeric(IDL), Cangrejo = as.numeric(Cangrejo)) Una vez se han leído y depurado estos dos conjuntos, se puede hacer la unión de ambos. Este nuevo conjunto es con el que se trata para los métodos de clasificación. occ_mof <- occ %> % inner_join(mof, by=c("Cangrejo","IDL","fecha")) El fichero original de EV contenía 18 variables, sin embargo solo siete de ellas contenían información útil. El formato que tenían las coordenadas latitud y longitud no era numérico sino una cadena de texto con puntos tanto en los miles como para la coma decimal. Se corrige el formato de la cadena de texto y se transforma en numérico. La variable hábitat se transforma en factor y se cambian los nombres de sus categorías para que sean más cortos y estén en español. Por último, aquellos registros que tuvieran un valor perdido fueron descartados, ya que la única variable con valores perdidos era la que registra el número de capturas, y tener ese valor es la finalidad de este conjunto. ev_bruto <- read_delim("data/Don_biom_red-swamp-crayfish_ev_20230111.csv", delim = ";",escape_double = FALSE, col_types = cols(eventDate = col_date(format = " %d/ %m/ %Y"), eventTime = col_time(format = " %H: %M: %Sz"), decimalLatitude = col_character()), locale = locale(encoding = "ISO-8859-1"), trim_ws = TRUE) ev <- ev_bruto %> % select(eventDate, localityID, eventTime,decimalLatitude, decimalLongitude, habitat, sampleSizeValue) %> % 100 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO A.1. IMPORTACIÓN mutate( decimalLatitude = as.numeric(str_remove(decimalLatitude,"\\."))/100, decimalLongitude = as.numeric(decimalLongitude)/100) %> % mutate(habitat = as.factor(habitat) %> % `levels<-`(c("Lagunas temporales","Marismas"))) %> % drop_na() El último conjunto es el que registra las variables ambientales, como ya se comentó este proviene de un excel, a diferencia del resto que provenian de ficheros csv. Además, estos datos no fueron publicados al mismo tiempo que los otros, y por tanto requirieron de un trabajo de depuración más meticuloso. La variable originalmente nombrada Profundidad de medida...12 tenía solo dos posibles valores, o bien “Superficie” o bien NA, en cuyo caso ninguna de las variables numéricas tomaba un valor, eran todos perdidos también. Por tanto, lo primero que se hizo fue filtrar que esta variable no tuviera un valor perdido. Después de esto se eliminaron aquellas variables que no tenían información útil para los análisis posteriores. En este caso el número de variables de utilidad era mayor, por lo que se realizó la selección a la inversa. En varios de los casos las variables tenían espacios o caracteres extraños, por lo que se renombraron para evitar tildes, ñ y espacios. Se creó una nueva variable que agrupase Mes y Año, para poder hacer visualizaciones de las variables ambientales a lo largo del tiempo. También se corrigió el tipo de algunas que no fue registrado correctamente. Finalmente tuvo que aplicarse la función unique() pues la variable Código nasa generaba gran cantidad de registros por cada localización y mes sin ninguna variación en las variables ambientales, pero al eliminarla ya no era necesario tener tantas repeticiones. Además de esto, tras eliminar los valores repetidos se advirtió que algunos de los datos tenían un comportamiento extraño. Para una misma fecha y localización había muchos registros que solo diferían en una única variable, y esta variable iba aumentando 1 unidad su valor en cada registro. Se asumió que se trataba de un error al introducir los datos, ya que esto producía valores imposibles en el pH entre otras cosas. Por tanto fueron eliminados todos los pares Fecha y localización con 3 o más registros diferentes. pres_abs_bruto <- read_excel("data/Doñana_Crayfish_Databases.xlsx", sheet = "Presence-Absence + FQvariables", col_types = c("numeric","numeric","numeric", "numeric","text","text","text", "numeric","text","text","text", "text","text","numeric","numeric", "numeric","numeric","numeric", "numeric","numeric","numeric", "numeric","numeric"), na = "N/D") pres_abs_depurado <- pres_abs_bruto %> % filter(!is.na(`Profundidad de medida...12`)) %> % select(-ID, -`Tipo nasa`, -`Código nasa`, -`Procambarus clarkii`, -`Profundidad de medida...12`, -`Profundidad de medida...13`) %> % APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO 101 A.2. ANÁLISIS DESCRIPTIVO rename(Ano = Año, Habitat = Hábitat, Sup_encharcada = `Superficie encharcada`, Prof_max = `Profundidad máxima`, Vegetacion = Vegetación, Temperatura = `Temperatura del agua`, O2_porc = `Oxígeno disuelto en porcentaje`, O2_abs = `Oxígeno disuelto en valor absoluto`, Clorofila = `Clorofila a`) %> % mutate(Fecha = Ano + (Mes-1)/12, Habitat = as.factor(Habitat), Sup_encharcada = as.factor(Sup_encharcada), Prof_max = as.numeric(Prof_max), Vegetacion = as.factor(Vegetacion)) %> % unique() pres_abs <- pres_abs_depurado %> % anti_join(pres_abs_depurado %> % group_by(Ano, Mes, IDL) %> % summarise(registros = n()) %> % filter(registros > 2), by=c("Ano","Mes","IDL")) La visualización de estos datos no se hizo con summary() como en los conjuntos anteriores, sino que se recurrió a la librería skimr para obtener algunas medidas de las variables numéricas y estas se reflejaron en una tabla. pa_skim <- pres_abs %> % select_if(is.numeric) %> % skim() pa_tab <- cbind(pa_skim$n_missing, pa_skim$complete_rate, pa_skim$numeric.mean, pa_skim$numeric.sd) rownames(pa_tab) <- pres_abs %> % select_if(is.numeric) %> % names() pa_tab %> % round(3) %> % kable(col.names = c("Valores perdidos"," % completos", "Media","Desviación típica"), caption = "\\label{PresAbsNum}Variables ambientales numéricas") Para la visualización de las variables categóricas sí que se usó summary() tras seleccionar solo esas tres variables. A.2. Análisis descriptivo Una vez realizada la importación y limpieza de los datos, lo siguiente fue el análisis descriptivo. En la mayoría de los casos no hubo que hacer modificaciones a los conjuntos para 102 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO A.2. ANÁLISIS DESCRIPTIVO crear las gráficas, y el código de estas no tiene especial interés por lo que no será incluido en este apéndice. Sin embargo, para el conjunto EV sí que se realizó una transformación para hacer la ponderación y la posterior agrupación por trimestres. La ponderación consistió en hacer la agrupación por localización, hallar la media de las capturas en cada una y luego crear una variable que almacenase el cociente entre las capturas reales registradas y la media habitual en esa localización. La información se guardó en un nuevo objeto ev_graf ya que inicialmente se utilizó para el apartado de Análisis gráfico. ev %> % select(eventDate, localityID, sampleSizeValue, habitat) %> % group_by(localityID) %> % mutate(media = mean(sampleSizeValue, na.rm=TRUE), capt = sampleSizeValue/media) -> ev_graf Para la agrupación por trimestres se partió del conjunto ev_graf y se creó una variable para almacenar el trimestre al que correspondía el registro. ev_graf %> % mutate(fecha_trim = floor(month(eventDate)/3) /4+ year(eventDate)) %> % group_by(fecha_trim, localityID, habitat) %> % summarise(sampleSizeValue = sum(sampleSizeValue, na.rm = TRUE)) %> % ungroup() %> % group_by(localityID) %> % mutate(media = mean(sampleSizeValue, na.rm=TRUE), capt = sampleSizeValue/media, estacion = fecha_trim % %1*4+1, est_nombre = as.factor(c("Inv","Pri","Ver","Oto")[estacion]) )-> ev_graf_trimestre Para las tablas que reflejaban la media y desviación típica de las capturas en trampas en función del hábitat y la estación, se aplicó primero la función skim() de la librería skimr. Posteriormente se seleccionó la medida de interés, habiendo filtrado previamente solo aquellos datos que se correspondían a sampleSizeValue. Una vez seleccionados los datos, se reordenaron para ser visualizados de forma más sencilla y gracias a la función kable() de kableExtra se creó la tabla ya vista anteriormente. ev_graf_skim <- ev_graf_trimestre %> % ungroup() %> % group_by(est_nombre, habitat) %> % skim() %> % filter(skim_variable == "sampleSizeValue") ev_graf_t1 <- ev_graf_skim %> % select(est_nombre, habitat, numeric.mean) %> % rename(media = numeric.mean) %> % mutate(media = round(media,2)) %> % pivot_wider(names_from = c("est_nombre"), values_from = c("media")) ev_graf_t2 <- ev_graf_skim %> % select(est_nombre, habitat, numeric.sd) %> % APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO 103 A.3. MODELOS CLASIFICACIÓN rename(desviacion_tipica = numeric.sd) %> % mutate(desviacion_tipica = round(desviacion_tipica,2)) %> % pivot_wider(names_from = c("est_nombre"), values_from = c("desviacion_tipica")) A.3. Modelos clasificación Para predecir las variables Sexo y Madurez se utilizaron como variables predictoras las variables sobre el tamaño de los cangrejos: largo del cuerpo, largo del cefalotorax, ancho del cefalotorax y peso. Se crearon modelos por el método de los K-vecinos más cercanos y regresión logística tanto usando solo el largo del cuerpo como con las demás. Para árboles de decisión solo se usó el largo del cuerpo y para bosques aleatorios se utilizaron todas. Las muestras de entrenamiento y test se generaron de forma estratificada respecto de la variable objetivo, y se mantuvieron las mismas para todos los métodos realizados. Hay dos subconjuntos por cada variable: el primero es el que no tiene ningún valor NA, ni en la variable objetivo ni en las variables predictoras. El segundo solo se necesita que sea conocida la variable objetivo, ya que el largo del cuerpo es conocido en todos los cangrejos registrados. Las muestras generadas utilizan un 70% de los datos para la muestra de entrenamiento y un 30% para muestra test. occ_mof_completo_sexo <- occ_mof %> % filter(!is.na(sex),!is.na(largo_cefalotorax), !is.na(ancho_cefalotorax),!is.na(peso)) occ_mof_comp_sexo_split <- occ_mof_completo_sexo %> % initial_split(prop=0.7,strata = sex) occ_mof_comp_sexo_ent <- occ_mof_comp_sexo_split %> % training() occ_mof_comp_sexo_test <- occ_mof_comp_sexo_split %> % testing() occ_mof_largo_sexo <- occ_mof %> % filter(!is.na(sex)) occ_mof_largo_sexo_split <- occ_mof_largo_sexo %> % initial_split(prop=0.7,strata = sex) occ_mof_largo_sexo_ent <- occ_mof_largo_sexo_split %> % training() occ_mof_largo_sexo_test <- occ_mof_largo_sexo_split %> % testing() occ_mof_completo_lfst <- occ_mof %> % filter(!is.na(lifestage),!is.na(largo_cefalotorax), !is.na(ancho_cefalotorax),!is.na(peso)) occ_mof_comp_lfst_split <- occ_mof_completo_lfst %> % initial_split(prop=0.7,strata = lifestage) 104 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO A.3. MODELOS CLASIFICACIÓN occ_mof_comp_lfst_ent <- occ_mof_comp_lfst_split %> % training() occ_mof_comp_lfst_test <- occ_mof_comp_lfst_split %> % testing() occ_mof_largo_lfst <- occ_mof %> % filter(!is.na(lifestage)) occ_mof_largo_lfst_split <- occ_mof_largo_lfst %> % initial_split(prop=0.7,strata = lifestage) occ_mof_largo_lfst_ent <- occ_mof_largo_lfst_split %> % training() occ_mof_largo_lfst_test <- occ_mof_largo_lfst_split %> % testing() A continuación se muestra el código utilizado para crear los modelos en cada caso. El de las gráficas o las matrices de confusión es omitido ya que es repetitivo y tiene menor interés. Para el modelo KNN con todas las variables predictoras se utiliza la librería caret, se aplica un centrado y escalado de las variables predictoras y se generan 11 valores para el número de vecinos, quedando en el modelo final el que mejor rendimiento tuviera en la muestra entrenamiento. ctrl <- trainControl(method="cv",classProbs=TRUE, summaryFunction = defaultSummary ) KNN_todas_sexo <- train(sex ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_sexo_ent, method = "knn", trControl = ctrl, preProcess =c("center","scale"), tuneLength=11 ) KNN_todas_lfst <- train(lifestage ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_lfst_ent, method = "knn", trControl = ctrl, preProcess =c("center","scale"), tuneLength=11) Lo mismo se aplica para los modelos solo con el largo del cuerpo como variable predictora. KNN_largo_sexo <- train(sex ~ largo_cuerpo, data = occ_mof_largo_sexo_ent, APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO 105 A.3. MODELOS CLASIFICACIÓN method = "knn", trControl = ctrl, preProcess =c("center","scale"), tuneLength=11) KNN_largo_lfst <- train(lifestage ~ largo_cuerpo, data = occ_mof_largo_lfst_ent, method = "knn", trControl = ctrl, preProcess =c("center","scale"), tuneLength=11) Los árboles de clasificación se generan con la librería C50 que con una única instrucción genera el árbol, al que luego puede aplicarse un plot() para visualizarlo tal como se mostró en el apartado correspondiente. arbol_largo_sexo <- C5.0(sex ~ largo_cuerpo, data = occ_mof_largo_sexo_ent) arbol_largo_lfst <- C5.0(lifestage ~ largo_cuerpo, data = occ_mof_largo_lfst_ent) Para el método de bosques aleatorios se utilizó la librería randomForest. Se indicó que el parámetros mtry fuera igual a tres, para así en cada nodo usar tres variables predictoras de las cuatro disponibles. De este modo hay variabilidad en los árboles por no usar siempre todas las variables pero permite optimizar en función de tres de ellas para obtener los mejores resultados. El número de árboles generados se fijó en 100. bosque_todas_sexo <- randomForest(sex ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_sexo_ent, importance = TRUE, replace = TRUE, ntree=100,mtry=3) bosque_todas_lfst <- randomForest(lifestage ~ largo_cuerpo + largo_cefalotorax + ancho_cefalotorax + peso, data = occ_mof_comp_lfst_ent, importance = TRUE, replace = TRUE, ntree=100,mtry=3) Los modelos generados de regresión logística se hicieron con la función glm(), y se hicieron para una única variable predictora y para las cuatro. 106 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO A.4. MODELOS REGRESIÓN logi_largo_sexo <- glm(sex ~ largo_cuerpo,family = binomial, data = occ_mof_largo_sexo_ent) logi_todas_sexo <- glm(sex ~ largo_cuerpo + ancho_cefalotorax + largo_cefalotorax + peso, family = binomial, data = occ_mof_comp_sexo_ent) logi_largo_lfst <- glm(lifestage ~ largo_cuerpo,family = binomial, data = occ_mof_largo_lfst_ent) logi_todas_lfst <- glm(lifestage ~ largo_cuerpo + ancho_cefalotorax + largo_cefalotorax + peso, family = binomial, data = occ_mof_comp_lfst_ent) A.4. Modelos regresión Al inicio del capítulo sobre métodos de regresión se generaron algunas gráficas y se representaron las matrices de correlación de las variables de tamaño, ya que era relevante para justificar las posteriores decisiones. La decisión tomada fue separar la información de mof en dos conjuntos en función de si los registros se tomaron antes o después del 1 de enero de 2010. El conjunto anterior a esta fecha se nombró mof_f1 y el posterior mof_f2, el código a continuación genera los dos conjuntos y muestra la gráfica de la matriz de correlación de cada uno. mof_f1 <- mof %> % filter(fecha < "2010-01-01") mof_f1 %> % select(largo_cuerpo, largo_cefalotorax, ancho_cefalotorax, peso) %> % cor(use="pairwise.complete.obs") %> % corrplot(method="ellipse") mof_f2 <- mof %> % filter(fecha > "2010-01-01") mof_f2 %> % select(largo_cuerpo, largo_cefalotorax, ancho_cefalotorax, peso) %> % cor(use="pairwise.complete.obs") %> % corrplot(method="ellipse") El primer modelo generado era de regresión lineal, queriendo predecir el largo del cefalotorax en función del largo del cuerpo, y tomando el conjunto de datos sin dividir. Primero se filtraron los casos completos, luego se creó la división en muestra de entrenamiento (80% de los casos) y muestra test (20% de los casos). El error cometido se calculó mediante validación cruzada sobre la muestra de entrenamiento. La muestra test se utilizó para visualizar las prediccciones en comparación con los valores reales de la variable. APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO 107 A.4. MODELOS REGRESIÓN mof_larcf <- mof %> % filter(!is.na(largo_cefalotorax)) mof_larcf_split <- initial_split(mof_larcf, prop=0.8, strata=largo_cefalotorax) reg_mof_larcf =lm(largo_cefalotorax ~ largo_cuerpo, data=training(mof_larcf_split)) VC_mof_1_lcf <- cv.lm(largo_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_larcf_split), seed = 1909,plotit = FALSE, printit = FALSE) Este modelo tenía como finalidad recalcar la necesidad de separar los dos conjuntos, ya que las estimaciones obtenidas no eran correctas. Una vez se visualizó esto, se repitió el proceso con cada uno de los subconjuntos antes generados. mof_larcf_f1 <- mof_f1 %> % filter(!is.na(largo_cefalotorax)) mof_larcf_f1_split <- initial_split(mof_larcf_f1, prop=0.8, strata=largo_cefalotorax) reg_mof_larcf_f1 =lm(largo_cefalotorax ~ largo_cuerpo, data=training(mof_larcf_f1_split)) VC_mof_larcf_f1 <- cv.lm(largo_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_larcf_f1_split), seed = 1909,plotit = FALSE, printit = FALSE) mof_larcf_f2 <- mof_f2 %> % filter(!is.na(largo_cefalotorax)) mof_larcf_f2_split <- initial_split(mof_larcf_f2, prop=0.8, strata=largo_cefalotorax) reg_mof_larcf_f2 =lm(largo_cefalotorax ~ largo_cuerpo, data=training(mof_larcf_f2_split)) VC_mof_larcf_f2 <- cv.lm(largo_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_larcf_f2_split), seed = 1909,plotit = FALSE, printit = FALSE) El error cometido por estos modelos era considerablemente menor que el del modelo sin división temporal de los datos. También las predicciones eran notablemente más acertadas que en el otro modelo. Para la variable ancho del cefalotorax ocurría lo mismo, por tanto se generaron únicamente los modelos con los datos separados, sin hacer primero el modelo incorrecto. 108 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO A.4. MODELOS REGRESIÓN mof_anccf_f1 <- mof_f1 %> % filter(!is.na(ancho_cefalotorax)) mof_anccf_f1_split <- initial_split(mof_anccf_f1, prop=0.8, strata=ancho_cefalotorax) reg_mof_anccf_f1 =lm(ancho_cefalotorax ~ largo_cuerpo, data=training(mof_anccf_f1_split)) VC_mof_anccf_f1 <- cv.lm(ancho_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_anccf_f1_split), seed = 1909,plotit = FALSE, printit = FALSE) mof_anccf_f2 <- mof_f2 %> % filter(!is.na(ancho_cefalotorax)) mof_anccf_f2_split <- initial_split(mof_anccf_f2, prop=0.8, strata=ancho_cefalotorax) reg_mof_anccf_f2 =lm(ancho_cefalotorax ~ largo_cuerpo, data=training(mof_anccf_f2_split)) VC_mof_anccf_f2 <- cv.lm(ancho_cefalotorax ~ largo_cuerpo, m=10, data = training(mof_anccf_f2_split), seed = 1909,plotit = FALSE, printit = FALSE) Para la variable peso en la gráfica inicial no se veía que tuviera un comportamiento demasiado distinto entre antes y después de 2010, de modo que se creó primero el modelo usando el conjunto total. mof_pes <- mof %> % filter(!is.na(peso)) mof_pes_split <- initial_split(mof_pes, prop=0.8,strata=peso) reg_mof_peso =lm(peso ~ largo_cuerpo, data=training(mof_pes_split)) VC_mof_1_peso <- cv.lm(peso ~ largo_cuerpo, m=10, data = training(mof_pes_split), seed = 1909,plotit = FALSE, printit = FALSE) Este modelo era mejor que en el caso del largo del cefalotorax usando todos los datos, pero era peor que aquellos modelos creados con los datos divididos. Igualmente el ajuste no daba la impresión de ser malo por el mismo motivo que ocurría con los tamaños del cefalotorax, pero por confirmarlo se crearon los modelos haciendo la separación de los datos. APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO 109 A.4. MODELOS REGRESIÓN knn_tune_model_reg <- nearest_neighbor(neighbors = tune(), weight_func = tune(), dist_power = tune()) %> % set_engine("kknn") %> % set_mode("regression") knn_tune_wkfl <- workflow() %> % add_model(knn_tune_model_reg) %> % add_formula(capt ~ .) knn_grid <- grid_random(parameters(knn_tune_model_reg), size = 10) lag_folds <- vfold_cv(training(lag_imp_split), v = 10, strata = capt) lag_knn_tuning <- knn_tune_wkfl %> % tune_grid(resamples = lag_folds, grid = knn_grid) lag_knn_best_model <- lag_knn_tuning %> % select_best(metric = "rmse") lag_knn_tune_wkfl_final <- knn_tune_wkfl %> % finalize_workflow(lag_knn_best_model) mar_folds <- vfold_cv(training(mar_imp_split), v = 10, strata = capt) mar_knn_tuning <- knn_tune_wkfl %> % tune_grid(resamples = mar_folds, grid = knn_grid) mar_knn_best_model <- mar_knn_tuning %> % select_best(metric = "rmse") mar_knn_tune_wkfl_final <- knn_tune_wkfl %> % finalize_workflow(mar_knn_best_model) Para las lagunas el valor óptimo de k era 12 y para las marismas 13. En ambos casos el error cometido fue similar al de sus bosques aleatorios equivalentes, por lo que se concluyó que era el mejor resultado que podía obtenerse con los datos disponibles. 116 APÉNDICE A. APÉNDICE: CÓDIGO OMITIDO Apéndice B Apéndice: Código descartado A lo largo de la realización de este trabajo, se plantearon múltiples ideas que no obtuvieron buenos resultados. Algunas de ellas quedan reflejadas en este apéndice, para justificar que no fueran llevadas a cabo en la parte principal del documento. B.1. Regresión para capturas en trampas B.1.1. Datos sin imputar, por hábitat Hay 24 registros en la muestra test y 11 predictores numéricos. Si se usan todas en un único modelo, el resultado obtenido sería: summary(lm(capt ~ ., data = lag_sin_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = capt ~ ., data = lag_sin_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max ## -0.32614 -0.03342 0.01071 0.02713 0.22858 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) -1.132e+02 1.183e+02 -0.957 0.357 ## Anyo 5.681e-02 5.883e-02 0.966 0.353 ## Mes -1.654e-03 1.376e-02 -0.120 0.906 ## Prof_max -1.744e-02 3.366e-02 -0.518 0.614 ## Temperatura 4.137e-02 3.794e-02 1.090 0.297 ## Conductividad -6.410e-02 5.215e-02 -1.229 0.243 ## O2_porc -1.177e-02 5.084e-02 -0.232 0.821 ## pH -1.413e-02 5.303e-02 -0.267 0.794 ## Turbidez -1.439e-03 6.148e-02 -0.023 0.982 117 B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS Tabla B.1: Correlación capturas lagunas respecto al clima sin imputar Conductividad Temperatura Nitrato pH Mes O2_porc 0.262 0.221 0.176 0.16 0.137 0.117 Tabla B.2: Regresiones capturas lagunas, clima sin imputar Variables R^2 ajustado p-valor Conductividad + Temperatura + Nitrato + pH + Mes + O2_porc -0.01725 0.4955 Conductividad + Temperatura + Nitrato + pH + Mes 0.03892 0.35454 Conductividad + Temperatura + Nitrato + pH 0.07191 0.25777 Conductividad + Temperatura + Nitrato 0.09353 0.1813 Conductividad + Temperatura 0.03879 0.25396 ## Clorofila -4.855e-03 5.702e-02 -0.085 0.934 ## Amonio -1.645e-02 8.469e-02 -0.194 0.849 ## Nitrato 4.658e-02 1.264e-01 0.368 0.719 ## ## Residual standard error: 0.1373 on 12 degrees of freedom ## Multiple R-squared: 0.3101, Adjusted R-squared: -0.3224 ## F-statistic: 0.4902 on 11 and 12 DF, p-value: 0.876 El p-valor es muy alto, y el R2ajustado llega a ser negativo. Hay demasiadas variables predictoras, e irlas quitando una a una sería excesivamente lento. Una opción es partir del modelo con aquellas variables que tengan mayor correlación con la variable objetivo, y luego irlo refinando. capt_lag_sin_cor <- lag_sin_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") capt_lag_sin_cor2 <- capt_lag_sin_cor["capt",] %> % abs() %> % sort(decreasing = TRUE) capt_lag_sin_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación capturas lagunas respecto al clima sin imputar") Partiendo de las variables con mayor correlación, se va reduciendo el número de variables del modelo y se muestra en la Tabla B.2 los resultados. En todos estos modelos se acepta que los coeficientes sean 0 para el modelo total, y el coeficiente de determinación apenas llega a 0.1. Se calcula el porcentaje de error en el modelo con mejor coeficiente de determinación: err_reg_ev_lag_sin <- lag_sin_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Conductividad + Temperatura + Nitrato, data = lag_sin_imp_ent), lag_sin_imp_test)) %> % 118 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS rmse(truth = capt, estimate = pred) %> % pluck(3) lag_sin_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Conductividad + Temperatura + Nitrato, data = lag_sin_imp_ent), lag_sin_imp_test)) %> % ggplot(aes(x=capt, y = pred)) + geom_point() + geom_abline() + theme_minimal() 0.9 1.0 1.1 1.2 0.9 1.0 1.1 1.2 capt pred ## [1] "El error cometido es del 17.43 %" Esta misma estructura se seguirá en las marismas, también en ambos hábitats para los datos del clima imputados. En las marismas el número de registros para la muestra entrenamiento es 30. summary(lm(capt ~ ., data = mar_sin_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = capt ~ ., data = mar_sin_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 119 B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS Tabla B.3: Correlación capturas marismas respecto al clima sin imputar Amonio Prof_max O2_porc Mes Anyo pH 0.449 0.289 0.223 0.198 0.169 0.149 ## -0.19855 -0.08196 0.01032 0.07676 0.19414 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) -62.359527 69.759444 -0.894 0.38316 ## Anyo 0.031500 0.034693 0.908 0.37590 ## Mes 0.015526 0.011053 1.405 0.17712 ## Prof_max -0.063786 0.036717 -1.737 0.09943 . ## Temperatura 0.026075 0.041001 0.636 0.53280 ## Conductividad -0.007796 0.028223 -0.276 0.78551 ## O2_porc 0.086723 0.038843 2.233 0.03851 * ## pH -0.007216 0.043264 -0.167 0.86939 ## Turbidez 0.009176 0.040316 0.228 0.82252 ## Clorofila -0.030165 0.039769 -0.758 0.45797 ## Amonio -0.103100 0.029998 -3.437 0.00294 ** ## Nitrato -0.052303 0.034547 -1.514 0.14740 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 ## ## Residual standard error: 0.1505 on 18 degrees of freedom ## Multiple R-squared: 0.5764, Adjusted R-squared: 0.3175 ## F-statistic: 2.227 on 11 and 18 DF, p-value: 0.06358 El modelo total tiene mejor ajuste que en las lagunas, y aunque a un 5 % de confianza siga sin rechazarse la nulidad de coeficientes, a un 10 % sí se rechazaría. capt_mar_sin_cor <- mar_sin_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") capt_mar_sin_cor2 <- capt_mar_sin_cor["capt",] %> % abs() %> % sort(decreasing = TRUE) capt_mar_sin_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación capturas marismas respecto al clima sin imputar") En este caso al menos se rechaza que los coeficientes del modelo sean nulos, aunque el coeficiente de determinación dista mucho de ser suficientemente grande. Se mide el RMSE de la muestra test en el modelo con mejor coeficiente de determinación. 120 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS Tabla B.4: Regresiones capturas marismas, clima sin imputar Variables R^2 ajustado p-valor Amonio + Prof_max + O2_porc + Mes + Anyo + pH 0.33997 0.01334 Amonio + Prof_max + O2_porc + Mes + Anyo 0.35558 0.00696 Amonio + Prof_max + O2_porc + Mes 0.35904 0.00397 Amonio + Prof_max + O2_porc 0.36206 0.00201 Amonio + Prof_max 0.21966 0.01339 err_reg_ev_mar_sin <- mar_sin_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Amonio + Prof_max + O2_porc, data = mar_sin_imp_ent), mar_sin_imp_test)) %> % rmse(truth = capt, estimate = pred) %> % pluck(3) mar_sin_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Amonio + Prof_max + O2_porc, data = mar_sin_imp_ent), mar_sin_imp_test)) %> % ggplot(aes(x=capt, y = pred)) + geom_point() + geom_abline() + theme_minimal() 0.9 1.0 1.1 1.2 1.0 1.1 1.2 capt pred ## [1] "El error cometido es del 15.96 %" APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 121 B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS B.1.2. Datos imputados, por hábitat summary(lm(capt ~ ., data = lag_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = capt ~ ., data = lag_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max ## -0.81339 -0.08796 0.00142 0.09074 1.24595 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) -4.266531 26.264679 -0.162 0.8713 ## Anyo 0.002642 0.013070 0.202 0.8402 ## Mes -0.007997 0.008416 -0.950 0.3443 ## Prof_max 0.012443 0.027814 0.447 0.6556 ## Temperatura 0.053084 0.035598 1.491 0.1390 ## Conductividad 0.007738 0.026787 0.289 0.7733 ## O2_abs 0.028750 0.031360 0.917 0.3614 ## pH -0.075226 0.027179 -2.768 0.0067 ** ## Turbidez -0.069947 0.046932 -1.490 0.1392 ## Clorofila 0.029030 0.026653 1.089 0.2787 ## Amonio 0.053832 0.034474 1.561 0.1215 ## Nitrato 0.004219 0.031265 0.135 0.8929 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 ## ## Residual standard error: 0.2721 on 102 degrees of freedom ## Multiple R-squared: 0.1162, Adjusted R-squared: 0.02085 ## F-statistic: 1.219 on 11 and 102 DF, p-value: 0.2843 capt_lag_imp_cor <- lag_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") capt_lag_imp_cor2 <- capt_lag_imp_cor["capt",] %> % abs() %> % sort(decreasing = TRUE) capt_lag_imp_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación capturas lagunas respecto al clima imputado") 122 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS Tabla B.5: Correlación capturas lagunas respecto al clima imputado pH Turbidez Mes Temperatura O2_abs Clorofila 0.205 0.113 0.091 0.078 0.068 0.058 Tabla B.6: Regresiones capturas lagunas, clima imputado Variables R^2 ajustado p-valor pH + Turbidez + Mes + Temperatura + O2_abs + Clorofila 0.04121 0.10399 pH + Turbidez + Mes + Temperatura + O2_abs 0.04415 0.07817 pH + Turbidez + Mes + Temperatura 0.03725 0.0866 pH + Turbidez + Mes 0.03692 0.06789 pH + Turbidez 0.0416 0.03511 err_reg_ev_lag_imp <- lag_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ pH + Turbidez, data = lag_imp_ent), lag_imp_test)) %> % rmse(truth = capt, estimate = pred) %> % pluck(3) lag_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ pH + Turbidez, data = lag_imp_ent), lag_imp_test)) %> % ggplot(aes(x=capt, y = pred)) + geom_point() + geom_abline() + theme_minimal() 0.6 0.7 0.8 0.9 1.0 1.1 0.5 1.0 1.5 capt pred APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 123 B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS ## [1] "El error cometido es del 27 %" summary(lm(capt ~ ., data = mar_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = capt ~ ., data = mar_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max ## -0.78389 -0.09273 0.01632 0.14610 0.86428 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) -54.202606 23.994444 -2.259 0.0259 * ## Anyo 0.027480 0.011939 2.302 0.0233 * ## Mes -0.003914 0.008006 -0.489 0.6259 ## Prof_max -0.040595 0.031627 -1.284 0.2021 ## Temperatura 0.048820 0.028982 1.684 0.0950 . ## Conductividad -0.011699 0.024619 -0.475 0.6356 ## O2_porc -0.008273 0.026213 -0.316 0.7529 ## pH 0.032059 0.029882 1.073 0.2857 ## Turbidez -0.055311 0.032258 -1.715 0.0893 . ## Clorofila 0.047156 0.030944 1.524 0.1305 ## Amonio -0.021830 0.027188 -0.803 0.4238 ## Nitrato 0.023537 0.030246 0.778 0.4382 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 ## ## Residual standard error: 0.2628 on 107 degrees of freedom ## Multiple R-squared: 0.1353, Adjusted R-squared: 0.04641 ## F-statistic: 1.522 on 11 and 107 DF, p-value: 0.134 capt_mar_imp_cor <- mar_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") capt_mar_imp_cor2 <- capt_mar_imp_cor["capt",] %> % abs() %> % sort(decreasing = TRUE) capt_mar_imp_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación capturas marismas respecto al clima imputado") err_reg_ev_mar_imp <- mar_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Turbidez + Anyo + Nitrato, 124 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.1. REGRESIÓN PARA CAPTURAS EN TRAMPAS Tabla B.7: Correlación capturas marismas respecto al clima imputado Turbidez Anyo Nitrato Clorofila pH Mes 0.189 0.185 0.124 0.099 0.079 0.073 Tabla B.8: Regresiones capturas marismas, clima imputado Variables R^2 ajustado p-valor Turbidez + Anyo + Nitrato + Clorofila + pH + Mes 0.04323 0.08886 Turbidez + Anyo + Nitrato + Clorofila + pH 0.05052 0.05357 Turbidez + Anyo + Nitrato + Clorofila 0.05224 0.03821 Turbidez + Anyo + Nitrato 0.05214 0.02725 Turbidez + Anyo 0.04322 0.0286 data = mar_imp_ent), mar_imp_test)) %> % rmse(truth = capt, estimate = pred) %> % pluck(3) mar_imp_test %> % select(capt) %> % cbind(pred = predict(lm(capt ~ Turbidez + Anyo + Nitrato, data = mar_imp_ent), mar_imp_test)) %> % ggplot(aes(x=capt, y = pred)) + geom_point() + geom_abline() + theme_minimal() 0.7 0.8 0.9 1.0 1.1 0.6 0.9 1.2 capt pred ## [1] "El error cometido es del 19.28 %" APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 125 B.2. PROPORCIÓN DE MACHOS Y HEMBRAS Tabla B.11: Correlación proporción marismas respecto al clima sin imputar Temperatura Mes Conductividad Prof_max O2_porc Amonio 0.375 0.214 0.133 0.11 0.074 0.051 Tabla B.12: Regresiones proporción marismas, clima sin imputar Variables R^2 ajustado p-valor Temperatura + Mes + Conductividad + Prof_max + O2_porc + Amonio 0.07437 0.34205 Temperatura + Mes + Conductividad + Prof_max + O2_porc -0.15243 0.77321 Temperatura + Mes + Conductividad + Prof_max -0.0782 0.63209 Temperatura + Mes + Conductividad -0.0121 0.45157 Temperatura + Mes 0.04234 0.26898 caption= "Correlación proporción marismas respecto al clima sin imputar") El pvalor se va reduciendo, pero incluso con un 10 % de significación, en todos los casos se acepta que los coeficientes sean nulos. err_reg_occ_mar_sin <- occ_mar_sin_imp_ent %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~Temperatura + Mes, data = occ_mar_sin_imp_ent))) %> % rmse(truth = proporcion, estimate = pred) %> % pluck(3) occ_mar_sin_imp_test %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~ Temperatura + Mes, data = occ_mar_sin_imp_ent), occ_mar_sin_imp_test)) %> % ggplot(aes(x=proporcion, y = pred)) + geom_point() + geom_abline() + theme_minimal() 132 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.2. PROPORCIÓN DE MACHOS Y HEMBRAS 45 50 55 30 40 50 60 70 proporcion pred ## [1] "El error cometido es del 23.77 %" B.2.2. Datos imputados, por hábitat summary(lm(proporcion ~ ., data = occ_lag_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = proporcion ~ ., data = occ_lag_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max ## -27.4777 -5.6715 0.1159 5.5986 28.8033 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 49.6798 1.4893 33.357 <2e-16 *** ## Ano 2.0684 2.4238 0.853 0.3976 ## Mes 0.7498 2.0352 0.368 0.7142 ## Prof_max -0.7697 1.9697 -0.391 0.6977 ## Temperatura 5.3478 2.4976 2.141 0.0373 * ## Conductividad -1.3794 1.3756 -1.003 0.3209 ## O2_abs 0.2038 1.7868 0.114 0.9096 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 133 B.2. PROPORCIÓN DE MACHOS Y HEMBRAS Tabla B.13: Correlación proporción lagunas respecto al clima imputado Temperatura Amonio Ano Clorofila Turbidez Prof_max 0.327 0.185 0.165 0.148 0.123 0.078 Tabla B.14: Regresiones proporción lagunas, clima imputado Variables R^2 ajustado p-valor Temperatura + Amonio + Ano + Clorofila + Turbidez + Prof_max 0.06363 0.14383 Temperatura + Amonio + Ano + Clorofila + Turbidez 0.07571 0.0956 Temperatura + Amonio + Ano + Clorofila 0.07801 0.07312 Temperatura + Amonio + Ano 0.0736 0.0617 Temperatura + Amonio 0.08325 0.03008 ## pH -1.0157 1.6985 -0.598 0.5526 ## Turbidez -3.1811 3.6010 -0.883 0.3813 ## Clorofila -1.4015 1.4344 -0.977 0.3333 ## Amonio -0.3144 2.8742 -0.109 0.9133 ## Nitrato 0.7716 1.9054 0.405 0.6873 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 ## ## Residual standard error: 10.68 on 49 degrees of freedom ## Multiple R-squared: 0.1807, Adjusted R-squared: -0.003208 ## F-statistic: 0.9826 on 11 and 49 DF, p-value: 0.4747 De nuevo no se rechaza la nulidad de coeficientes, y el R2ajustado es negativo. prop_lag_imp_cor <- occ_lag_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") prop_lag_imp_cor2 <- prop_lag_imp_cor["proporcion",] %> % abs() %> % sort(decreasing = TRUE) prop_lag_imp_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación proporción lagunas respecto al clima imputado") err_reg_occ_lag_imp <- occ_lag_imp_ent %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~Temperatura + Amonio, data = occ_lag_imp_ent))) %> % rmse(truth = proporcion, estimate = pred) %> % pluck(3) occ_lag_imp_test %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~ Temperatura + Amonio, data = occ_lag_imp_ent), 134 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.2. PROPORCIÓN DE MACHOS Y HEMBRAS occ_lag_imp_test)) %> % ggplot(aes(x=proporcion, y = pred)) + geom_point() + geom_abline() + theme_minimal() 42 46 50 54 40 50 60 70 proporcion pred ## [1] "El error cometido es del 19.9 %" summary(lm(proporcion ~ ., data = occ_mar_imp_ent %> % select_if(is.numeric))) ## ## Call: ## lm(formula = proporcion ~ ., data = occ_mar_imp_ent %>% select_if(is.numeric)) ## ## Residuals: ## Min 1Q Median 3Q Max ## -38.241 -7.090 2.877 7.836 37.359 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 49.12322 2.19108 22.420 <2e-16 *** ## Ano 2.88403 2.29036 1.259 0.214 ## Mes -2.08175 2.13681 -0.974 0.335 ## Prof_max -0.90283 2.08567 -0.433 0.667 ## Temperatura 1.55203 2.97622 0.521 0.604 ## Conductividad -0.69122 2.60381 -0.265 0.792 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 135 B.2. PROPORCIÓN DE MACHOS Y HEMBRAS Tabla B.15: Correlación proporción marismas respecto al clima imputado Nitrato Ano Clorofila Mes pH Amonio 0.203 0.173 0.114 0.113 0.062 0.055 Tabla B.16: Regresiones proporción marismas, clima imputado Variables R^2 ajustado p-valor Nitrato + Ano + Clorofila + Mes + pH + Amonio -0.00974 0.49934 Nitrato + Ano + Clorofila + Mes + pH 0.00791 0.37022 Nitrato + Ano + Clorofila + Mes 0.02454 0.25348 Nitrato + Ano + Clorofila 0.02337 0.23003 Nitrato + Ano 0.02357 0.18731 ## O2_abs 2.55248 2.74880 0.929 0.358 ## pH -0.09848 2.61046 -0.038 0.970 ## Turbidez -1.72825 4.45794 -0.388 0.700 ## Clorofila -2.00667 2.94319 -0.682 0.499 ## Amonio -0.37035 2.65271 -0.140 0.890 ## Nitrato -4.75425 2.96670 -1.603 0.115 ## --- ## Signif. codes: 0 ’***’ 0.001 ’**’ 0.01 ’*’ 0.05 ’.’ 0.1 ’ ’ 1 ## ## Residual standard error: 15.03 on 49 degrees of freedom ## Multiple R-squared: 0.1162, Adjusted R-squared: -0.08216 ## F-statistic: 0.5859 on 11 and 49 DF, p-value: 0.831 Igual que en las lagunas. prop_mar_imp_cor <- occ_mar_imp_ent %> % select_if(is.numeric) %> % cor(use = "pairwise.complete.obs") prop_mar_imp_cor2 <- prop_mar_imp_cor["proporcion",] %> % abs() %> % sort(decreasing = TRUE) prop_mar_imp_cor2[2:7] %> % t() %> % kable(digits=3, caption= "Correlación proporción marismas respecto al clima imputado") err_reg_occ_mar_imp <- occ_mar_imp_ent %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~ Nitrato + Ano, data = occ_mar_imp_ent))) %> % rmse(truth = proporcion, estimate = pred) %> % pluck(3) occ_mar_imp_test %> % select(proporcion) %> % cbind(pred = predict(lm(proporcion ~ Nitrato + Ano, 136 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.3. TUNING BOSQUES ALEATORIOS data = occ_mar_imp_ent), occ_mar_imp_test)) %> % ggplot(aes(x=proporcion, y = pred)) + geom_point() + geom_abline() + theme_minimal() 35 40 45 50 40 50 60 70 proporcion pred ## [1] "El error cometido es del 26.04 %" B.3. Tuning bosques aleatorios Inicialmente el tuning de árboles aleatorios no se llevaba a cabo por completo mediante tidymodels pues al indicar que el parámetro mtry se obtuviera mediante tuning el código producía un error. La forma de solucionarlo finalmente fue indicar los límites en los que debía buscarse el parámetro, pero antes de eso se aplicaba primero el tuning de mtry mediante caret y posteriormente utilizando el valor óptimo obtenido para mtry se aplicaba con tidymodels el tuning de trees y min_n. El código con el que se llevaba a cabo esto es el siguiente: lag_imp_rf_mtry <- train(capt ~ ., data=lag_imp_ent, method = 'rf', na.action = na.exclude, tuneGrid = expand.grid(.mtry = 3:10), trControl = control <- trainControl(method='repeatedcv', number=10, search = 'grid')) APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 137 B.3. TUNING BOSQUES ALEATORIOS ggplot(lag_imp_rf_mtry) + theme_minimal() + xlab("Número de variables seleccionadas") + ylab("RMSE") + ggtitle("Lagunas") rf_lag_capt_tune_model <- rand_forest(mtry = as.numeric(lag_imp_rf_mtry$bestTune), trees = tune(),min_n = tune()) %> % set_engine("randomForest") %> % set_mode("regression") lag_capt_rf_tune_wkfl <- workflow() %> % add_model(rf_lag_capt_tune_model) %> % add_formula(capt ~ .) lag_capt_folds <- vfold_cv(lag_imp_ent, v=10,strata = capt) rf_capt_grid <- grid_random(parameters(rf_lag_capt_tune_model),size=10) lag_capt_rf_tuning <- lag_capt_rf_tune_wkfl %> % tune_grid(resamples= lag_capt_folds, grid=rf_capt_grid, metrics = metric_set(rmse)) lag_capt_rf_best_model <- lag_capt_rf_tuning %> % select_best(metric = "rmse") lag_capt_rf_tune_wkfl_final <- lag_capt_rf_tune_wkfl %> % finalize_workflow(lag_capt_rf_best_model) lag_capt_rf_tune_err <- lag_capt_rf_tune_wkfl_final %> % last_fit(split = lag_imp_split) %> % collect_metrics() %> % cbind(media = mean(lag_imp_test$capt)) %> % summarise(error = .estimate/media*100) print(paste("El error cometido en las lagunas es del", round(lag_capt_rf_tune_err[1,1],2), " %")) 138 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO B.3. TUNING BOSQUES ALEATORIOS ## R version 4.2.1 (2022-06-23 ucrt) ## Platform: x86_64-w64-mingw32/x64 (64-bit) ## Running under: Windows 10 x64 (build 19045) ## ## Matrix products: default ## ## Random number generation: ## RNG: Mersenne-Twister ## Normal: Inversion ## Sample: Rounding ## ## locale: ## [1] LC_COLLATE=Spanish_Spain.utf8 LC_CTYPE=Spanish_Spain.utf8 ## [3] LC_MONETARY=Spanish_Spain.utf8 LC_NUMERIC=C ## [5] LC_TIME=Spanish_Spain.utf8 ## ## attached base packages: ## [1] grid stats graphics grDevices utils datasets methods ## [8] base ## ## other attached packages: ## [1] kknn_1.3.1 corrplot_0.92 DAAG_1.25.4 ## [4] magrittr_2.0.3 randomForest_4.7-1.1 C50_0.1.8 ## [7] rpart_4.1.16 caret_6.0-93 lattice_0.20-45 ## [10] yardstick_1.1.0 workflowsets_1.0.1 workflows_1.1.3 ## [13] tune_1.1.1 rsample_1.1.1 recipes_1.0.5 ## [16] parsnip_1.1.0 modeldata_1.1.0 infer_1.0.4 ## [19] dials_1.2.0 scales_1.2.1 broom_1.0.1 ## [22] tidymodels_1.0.0 vcd_1.4-11 kableExtra_1.3.4 ## [25] lubridate_1.8.0 knitr_1.42 readxl_1.4.1 ## [28] skimr_2.1.5 forcats_0.5.2 stringr_1.5.0 ## [31] dplyr_1.1.1 purrr_1.0.1 tidyr_1.3.0 ## [34] tibble_3.2.1 ggplot2_3.4.0 tidyverse_1.3.2 ## [37] readr_2.1.3 ## ## loaded via a namespace (and not attached): ## [1] backports_1.4.1 systemfonts_1.0.4 igraph_1.4.2 ## [4] plyr_1.8.7 repr_1.1.6 splines_4.2.1 ## [7] listenv_0.9.0 digest_0.6.29 foreach_1.5.2 ## [10] htmltools_0.5.4 fansi_1.0.3 googlesheets4_1.1.0 ## [13] tzdb_0.3.0 globals_0.16.2 modelr_0.1.11 ## [16] gower_1.0.0 vroom_1.6.1 svglite_2.1.1 ## [19] hardhat_1.3.0 jpeg_0.1-10 colorspace_2.0-3 ## [22] rvest_1.0.3 rbibutils_2.2.13 haven_2.5.1 ## [25] xfun_0.38 crayon_1.5.2 jsonlite_1.8.2 ## [28] libcoin_1.0-9 survival_3.3-1 zoo_1.8-12 ## [31] iterators_1.0.14 glue_1.6.2 gtable_0.3.3 ## [34] gargle_1.4.0 ipred_0.9-13 webshot_0.5.4 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO 139 B.3. TUNING BOSQUES ALEATORIOS ## [37] future.apply_1.10.0 mvtnorm_1.1-3 DBI_1.1.3 ## [40] Rcpp_1.0.9 viridisLite_0.4.1 Cubist_0.4.2.1 ## [43] GPfit_1.0-8 bit_4.0.5 proxy_0.4-27 ## [46] Formula_1.2-5 stats4_4.2.1 lava_1.7.2.1 ## [49] prodlim_2019.11.13 httr_1.4.5 RColorBrewer_1.1-3 ## [52] ellipsis_0.3.2 pkgconfig_2.0.3 farver_2.1.1 ## [55] deldir_1.0-6 nnet_7.3-17 dbplyr_2.3.2 ## [58] utf8_1.2.2 tidyselect_1.2.0 labeling_0.4.2 ## [61] rlang_1.1.0 DiceDesign_1.9 reshape2_1.4.4 ## [64] munsell_0.5.0 cellranger_1.1.0 tools_4.2.1 ## [67] cli_3.6.1 generics_0.1.3 evaluate_0.20 ## [70] fastmap_1.1.0 yaml_2.3.5 ModelMetrics_1.2.2.2 ## [73] bit64_4.0.5 fs_1.5.2 future_1.32.0 ## [76] nlme_3.1-157 xml2_1.3.3 compiler_4.2.1 ## [79] rstudioapi_0.14 png_0.1-8 e1071_1.7-11 ## [82] reprex_2.0.2 lhs_1.1.6 stringi_1.7.8 ## [85] Matrix_1.5-3 vctrs_0.6.1 pillar_1.9.0 ## [88] lifecycle_1.0.3 furrr_0.3.1 Rdpack_2.4 ## [91] lmtest_0.9-40 data.table_1.14.2 latticeExtra_0.6-30 ## [94] R6_2.5.1 parallelly_1.35.0 codetools_0.2-19 ## [97] MASS_7.3-58.1 withr_2.5.0 mgcv_1.8-42 ## [100] parallel_4.2.1 hms_1.1.3 modelenv_0.1.1 ## [103] timeDate_4022.108 class_7.3-20 rmarkdown_2.21 ## [106] inum_1.0-5 googledrive_2.1.0 partykit_1.2-20 ## [109] pROC_1.18.0 base64enc_0.1-3 interp_1.1-4 140 APÉNDICE B. APÉNDICE: CÓDIGO DESCARTADO Bibliografía [1] Alcorlo,Paloma yDiéguez-Uribeondo,Javier (2014). «El cangrejo señal y el declive de las poblaciones de cangrejo autóctono». Ambienta: La revista del Ministerio de Medio Ambiente,(109), pp. 52–61. [2] Alcorlo,Paloma;Geiger,Walter yOtero,Marina (2008). «Reproductive biology and life cycle of the invasive crayfish Procambarus clarkii (Crustacea: Decapoda) in diverse aquatic habitats of South-Western Spain: Implications for population control». Fundamental and Applied Limnology,173(3), p. 197. [3] Alcorlo,Paloma;Noguerales,Victor;Mollá,Salvador et al. (2019). «Modeling the population dynamics of the Red Swamp Crayfish (Procambarus clarkii) in Doñana: Application to the harvesting strategies». Journal of FisheriesSciences. com,13(2), pp. 001–011. [4] Blanco,Juan A. (2008). «Pasos básicos en la modelización ecológica», 33(Especial 2014), pp. 471–484. ISSN 2078-7235. [5] Bravo,Miguel A;Andreu,Ana C;Román,Isidro;Arribas,Rosa; Márquez-Ferrando,Rocío;Díaz-Delgado,Ricardo yBustamante,Javier (2023). «Long-term monitoring in the biometrics of the red-swamp crayfish (Procambarus clarkii, Girard 1852) in Doñana Wetlands 2004-2022». [6] Colomer,MªÀngels;Margalida,Antoni;Valencia,Luis yPalau,Anto- ni (2014). «Application of a computational model for complex fluvial ecosystems: The population dynamics of zebra mussel Dreissena polymorpha as a case study». Ecological Complexity,20, pp. 116–126. [7] Desconocido (2000). «Linear regression». https://web.archive.org/web/20080222195200/http://www.curvefit.com/linear_ regression.htm. [8] —— (2018). «Conejos, Cangrejos y Camalote». http://admreportaje.blogspot.com/2018/02/conejos-cangrejos-y-camalote.html. [9] —— (2023). «Random Forest: Bosque aleatorio. Definición y funcionamiento». https://datascientest.com/es/random-forest-bosque-aleatorio-definicion-y-funcionamiento. [10] Elin Waring,Amelia McNamara, Michael Quinn yZhu,Hao (2022). skimr: Compact and Flexible Summaries of Data. https://cran.r-project.org/web/packages/skimr/. R package version 2.1.5. 141