Revista Athenea · Artículo de investigaciónVol. 7 · Núm. 26 · Octubre–diciembre de 2026 · pp. 50–59 · ISSN-e 2737-6419
Página 50. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)Artículo de investigaciónhttps://doi.org/10.47460/athenea.v7i26.157Modelización no estacionaria de precipitaciones extremas medianteteoría de valores extremos y datos meteorológicosMaría de Lourdes Palacios*https://orcid.org/0000-0002-6149-2948mapalacios@espoch.edu.ecEscuela Superior Politécnica de ChimborazoRiobamba, EcuadorMartha Ximena Dávalos Villegashttps://orcid.org/0000-0001-7865-6307martha.davalos@espoch.edu.ecEscuela Superior Politécnica de ChimborazoRiobamba, EcuadorGuido Mazón-Fierrohttps://orcid.org/0000-0001-8745-2373gmazon@espoch.edu.ecEscuela Superior Politécnica de ChimborazoRiobamba, EcuadorRodrigo Toasahttps://orcid.org/0000-0001-8744-0794rodrigo.toasa@unach.edu.ecUniversidad Nacional de ChimborazoRiobamba, Ecuador*Autor de correspondencia: mapalacios@espoch.edu.ecRecibido: (13/07/2026), Aceptado: (20/10/2026)Resumen. La estimación de precipitaciones extremas suele asumir estacionariedad, aunque la variabil-idad climática puede modificar su distribución temporal. Este estudio evaluó la no estacionariedad demáximos anuales de precipitación diaria en tres regiones climáticas del Ecuador mediante modelos devalores extremos. Se procesaron series horarias de ERA5-Land, se incorporó el Índice Oceánico NiñoRelativo como covariable y se compararon formulaciones estacionarias y no estacionarias mediante máxi-ma verosimilitud, criterios de información, pruebas de razón de verosimilitudes y bootstrap paramétrico.Guayaquil mostró dependencia conjunta del tiempo y de la variabilidad asociada a El Niño-Oscilacióndel Sur; Puyo presentó una señal temporal sensible al producto de precipitación; Quito permaneciócompatible con estacionariedad. El control de robustez con CHIRPS mostró que la evidencia de cambiono es espacialmente uniforme y debe interpretarse considerando la incertidumbre de los productos enmalla.Palabras clave: valores extremos, precipitación, no estacionariedad, variabilidad climática.Non-stationary modeling of extreme precipitation using extreme value theory andmeteorological dataAbstract. Extreme precipitation estimation commonly assumes stationarity, although climate variabili-ty can modify its temporal distribution. This study assessed nonstationarity in annual maxima of dailyprecipitation across three climatic regions of Ecuador using extreme-value models. Hourly ERA5-Landseries were processed, the Relative Oceanic Niño Index was incorporated as a covariate, and stationa-ry and nonstationary formulations were compared through maximum likelihood, information criteria,likelihood-ratio tests, and parametric bootstrap. Guayaquil showed joint dependence on time and ElNiño-Southern Oscillation variability; Puyo exhibited a temporal signal sensitive to the precipitationproduct; Quito remained compatible with stationarity. A CHIRPS robustness analysis showed that evi-dence of change is not spatially uniform and should be interpreted together with uncertainty associatedwith gridded precipitation products.Keywords: extreme values, precipitation, nonstationarity, climate variability.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...50
Página 51. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)I.INTRODUCCIÓNLa precipitación extrema constituye una variable crítica para el diseño hidráulico, la gestión delriesgo, la planificación urbana y la evaluación de amenazas hidrometeorológicas. En estas aplicaciones,la teoría de valores extremos permite modelar probabilísticamente observaciones raras y estimar cuantilesasociados a eventos de baja frecuencia. En hidrología, los máximos por bloques y la distribución gen-eralizada de valores extremos constituyen herramientas consolidadas para representar máximos anualesde precipitación y caudal [1], [2].La formulación clásica supone que los parámetros de la distribución permanecen invariantes enel tiempo. Sin embargo, la variabilidad climática, los cambios en la circulación atmosférica y otrosprocesos pueden introducir tendencias o dependencias con covariables que contradigan ese supuesto.La discusión sobre la pérdida de estacionariedad en series hidrológicas ha impulsado modelos cuyosparámetros dependen del tiempo o de forzantes climáticos [3], [4]. A la vez, se ha advertido que lacomplejidad adicional puede amplificar la incertidumbre si la estructura no estacionaria se seleccionaúnicamente por mejor ajuste, sin evidencia física o comparaciones de parsimonia [5].En Ecuador, la evaluación de extremos de precipitación enfrenta además una marcada heterogenei-dad regional. La influencia del océano Pacífico es especialmente relevante en la Costa, mientras quela orografía andina y el transporte de humedad amazónica generan regímenes distintos en la Sierra yAmazonía. Estudios recientes han mostrado que los eventos extremos de El Niño modifican fuertementelas anomalías de precipitación en zonas costeras y que la señal se atenúa o cambia con la elevación y ladistancia al litoral [6]. En los Andes ecuatorianos, la alta variabilidad espacio-temporal y el efecto de latopografía dificultan generalizar el comportamiento de los eventos intensos [7]. Más recientemente, seha documentado una evolución multidecadal de la relación entre ENSO y la precipitación del Ecuador,con cambios en los patrones asociados a los sectores oriental y central del Pacífico [8].Los antecedentes disponibles para Ecuador caracterizan anomalías de precipitación durante even-tos El Niño, patrones espacio-temporales de extremos andinos y cambios multidecadales en la relaciónENSO-precipitación Paliz [6], [7], [8]. Sin embargo, no resuelven explícitamente si abandonar la esta-cionariedad mejora la representación GEV de máximos diarios ni cuánto modifica los niveles de retornode forma consistente entre regiones climáticamente contrastantes. Además, detectar una tendencia enuna serie de máximos no implica por sí solo que un modelo probabilístico no estacionario sea preferible.Por ello, resulta necesario comparar formalmente estructuras estacionarias y no estacionarias, cuantificarsu incertidumbre y evaluar la sensibilidad de las conclusiones al producto de precipitación empleado.El objetivo de este estudio fue evaluar la no estacionariedad de las precipitaciones extremas mediantemodelos de valores extremos que incorporan variación temporal y climática, comparando su desempeñocon una formulación estacionaria y cuantificando sus efectos sobre niveles de retorno en Guayaquil, Quitoy Puyo. La contribución no consiste en proponer una nueva distribución, sino en evaluar de manerahomogénea y parsimoniosa el supuesto de estacionariedad mediante una serie larga de reanálisis, unacovariable ENSO actualizada, selección multimodelo, bootstrap paramétrico y un control explícito derobustez con CHIRPS.El artículo se organiza de la siguiente manera. La Sección II resume la formulación de valoresextremos, la no estacionariedad y los criterios de selección. La Sección III describe los datos, el pre-procesamiento y el procedimiento estadístico. La Sección IV presenta y discute los resultados para lastres localidades, incluidos los niveles de retorno y el control de robustez. Finalmente, se exponen lasconclusiones y las implicaciones del análisis.II.FUNDAMENTOS TEÓRICOSA. Teoría de valores extremos y no estacionariedadPara una serie de máximos anuales, la distribución generalizada de valores extremos (GEV) integraen una única familia los tres dominios límite clásicos. Su función de distribución acumulada puedeexpresarse en términos de los parámetros de localización µ, escala σ > 0 y forma ξ [2]:G(z) = exp(−1 + ξz −µσ−1/ξ)(1)Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...51
Página 52. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)La condición de soporte exige 1 + ξ(z −µ)/σ > 0. En el caso estacionario, los tres parámetros per-manecen constantes. En una formulación no estacionaria, al menos uno puede vincularse a covariables.Para evitar sobreparametrización en series de longitud moderada, se mantuvieron constantes σ y ξ y sepermitió variación únicamente en µ, siguiendo una estructura lineal general:µt = β0 + β1Tt + β2Rt(2)donde Tt representa el tiempo estandarizado y Rt el índice climático estandarizado. Esta formulacióncontiene como casos particulares M0 (β1 = β2 = 0), M1 (β2 = 0), M2 (β1 = 0) y M3 (β1 y β2libres). El enfoque con covariables es consistente con la extensión de la teoría de extremos a problemashidrológicos no estacionarios [1], [4], pero requiere controlar la complejidad y la incertidumbre delmodelo [5].B. Selección, incertidumbre y niveles de retornoLos modelos se compararon mediante el criterio de información de Akaike corregido para muestrasfinitas (AICc), complementado con el criterio bayesiano de información (BIC) y pruebas de razónde verosimilitudes para estructuras anidadas.El AICc penaliza la complejidad adicional con mayorintensidad cuando el tamaño de muestra es limitado [9] y ha mostrado utilidad en la selección demodelos GEV estacionarios y no estacionarios de precipitación extrema [10]:AICc = 2k −2ℓ+ 2k(k + 1)n −k −1(3)La selección no se basó de manera automática en el menor AICc. Se consideraron también ∆AICc,BIC, significancia de modelos anidados, coherencia de los coeficientes y parsimonia. La incertidumbrede parámetros y cuantiles se estimó mediante bootstrap paramétrico, una estrategia adecuada cuandose requiere propagar la incertidumbre del ajuste hacia cantidades derivadas [11].Para una probabilidad anual de excedencia 1/T, el nivel de retorno estacionario asociado al periodoT se obtiene como:zT = µ + σξ(−ln1 −1T−ξ−1)(4)Cuando µ varía con el tiempo o con una covariable, zT deja de ser un cuantil fijo y debe interpretarsecondicionado al estado de las covariables. En consecuencia, los niveles calculados para los modelosno estacionarios se reportaron para años y valores específicos de RONI, evitando asumir que un únicoperiodo de retorno resume un proceso cambiante [12].III.METODOLOGÍAA. Datos y preprocesamientoSe seleccionaron tres localidades representativas de regímenes contrastantes del Ecuador continen-tal: Guayaquil en la Costa, Quito en los Andes y Puyo en la Amazonía. La selección respondió a trescriterios: representar los principales dominios climáticos del Ecuador continental, incorporar regímenesde precipitación controlados por mecanismos atmosféricos contrastantes y mantener un número reducidode localidades que permitiera aplicar un procedimiento homogéneo sobre series largas. Guayaquil repre-senta la influencia costera del Pacífico y de ENSO, Quito un régimen andino condicionado por la altitudy Puyo un ambiente húmedo amazónico. La Fig. 1 muestra su distribución espacial, mientras que laTabla 1 resume las características de las series empleadas.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...52
Página 53. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)Fig. 1. Localización de las áreas de estudio en las regiones Costa, Andes y Amazonía delEcuador.La Tabla 1 sintetiza la ubicación, el régimen climático representado, la fuente de información,el periodo analizado y las características de resolución y cobertura de las series utilizadas. Las treslocalidades fueron procesadas con el mismo producto, periodo temporal y procedimiento de extracción,lo que mantiene la comparabilidad del análisis.Tabla 1. Localidades de estudio y características de las series de precipitación.LocalidadRegiónCoordenadasFuente y periodoResolución ycoberturaGuayaquilCosta2,1710◦S; 79,9224◦OERA5-Land, 1951–20250,1◦; horaria; 27 394días; 75 máximosQuitoAndes0,1807◦S; 78,4678◦OERA5-Land, 1951–20250,1◦; horaria; 27 394días; 75 máximosPuyoAmazonía1,4924◦S; 77,9960◦OERA5-Land, 1951–20250,1◦; horaria; 27 394días; 75 máximosLa precipitación horaria se obtuvo de ERA5-Land mediante la colección ECMWF/ERA5_LAND/HOURLY y la banda total_precipitation_hourly, con resolución de 0,1◦y disponibilidad horaria desde1950 [13]. La extracción se ejecutó en Google Earth Engine [14]. Para cada localidad se utilizó elpíxel correspondiente a las coordenadas de referencia. Los acumulados diarios se construyeron en horalocal del Ecuador continental (UTC-5), integrando desde las 05:00 UTC hasta las 04:59 UTC del díasiguiente, con el fin de evitar la fragmentación de eventos nocturnos entre dos días UTC. Los valores setransformaron a milímetros; tres valores negativos de magnitud inferior a 10−5 mm se trataron comocero por corresponder a artefactos numéricos. La continuidad se verificó mediante el conteo de fechasúnicas por año, considerando 365 o 366 días según correspondiera, y mediante la búsqueda de registrosnulos y duplicados. No se identificaron días faltantes en las tres localidades; por tanto, se extrajo elmáximo diario de cada año calendario, obteniéndose 75 observaciones por localidad.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...53
Página 54. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)Como covariable climática se utilizó el Índice Oceánico Niño Relativo (RONI) del Climate PredictionCenter de NOAA, disponible desde 1950. RONI es una media móvil trimestral del índice relativo Niño 3.4y fue adoptado por el CPC en 2026 para el seguimiento operativo de ENSO [15]. A cada máximo anualse asignó el trimestre móvil aproximadamente centrado en el mes del evento (enero=DJF, febrero=JFM,. . . , diciembre=NDJ); luego, tiempo y RONI se estandarizaron por localidad para mejorar la estabilidadnumérica del ajuste. La asignación según el mes de ocurrencia permitió asociar cada máximo con elestado ENSO dominante alrededor del evento, evitando imponer un trimestre fijo a localidades condistinta estacionalidad pluviométrica.B. Análisis estadístico y modelizaciónLa exploración temporal se realizó con la prueba no paramétrica de Mann-Kendall [16] y la pen-diente de Sen [17]. Posteriormente se ajustaron por máxima verosimilitud cuatro modelos GEV: M0estacionario; M1 con µ dependiente del tiempo; M2 con µ dependiente de RONI; y M3 con µ depen-diente simultáneamente de tiempo y RONI. La escala y la forma permanecieron constantes en todaslas formulaciones. La comparación consideró log-verosimilitud, AICc, BIC, ∆AICc, pesos de Akaikey pruebas de razón de verosimilitudes. Los parámetros se estimaron minimizando la log-verosimilitudnegativa con la restricción σ > 0 y soporte GEV válido. Se utilizó el algoritmo L-BFGS-B y, cuandola convergencia requirió comprobación adicional, se contrastó la solución mediante Powell; se retuvola solución finita con mejor verosimilitud que cumpliera las restricciones del modelo. La incertidumbrese evaluó mediante 2 000 réplicas de bootstrap paramétrico por modelo seleccionado: en cada réplicase simuló una nueva muestra desde la GEV ajustada y se reajustó la misma estructura. Los intervalosdel 95 % correspondieron a los percentiles 2,5 y 97,5 de las distribuciones bootstrap. Los cálculos serealizaron en Python 3.13 con SciPy 1.17 y pandas 2.2, utilizando una semilla reproducible de 20261021.C. Control de robustez con CHIRPSSe añadió un análisis de sensibilidad con CHIRPS versión 3 [18] durante el periodo común 1998–2025. El intervalo 1998–2025 corresponde al periodo común con información disponible en ERA5-Landy en ambas variantes diarias de CHIRPS utilizadas; 1998 es el primer año disponible de DAILY_SATy 2025 el último año completo considerado.Se compararon las colecciones UCSB-CHC/CHIRPS/V3/DAILY_RNL y UCSB-CHC/CHIRPS/V3/DAILY_SAT. DAILY_RNL usa ERA5 y DAILY_SAT usaIMERG para distribuir temporalmente los mismos totales pentadales; por construcción, ninguna cons-tituye una validación diaria independiente y la documentación oficial recomienda contrastar ambascuando se analizan extremos diarios. CHIRPS conserva una convención de día UTC y no se desplazóa hora local, por lo que la comparación se interpretó como sensibilidad de los máximos anuales y nocomo correspondencia evento a evento. Se calcularon correlaciones de Pearson y Spearman, sesgo,error absoluto medio, raíz del error cuadrático medio, razón entre medias, concordancia de Willmott ypendientes de Sen.IV.RESULTADOS Y DISCUSIÓNA. Comportamiento temporal de los máximos anualesLos máximos anuales de precipitación diaria presentaron diferencias marcadas entre localidades.Guayaquil registró una media de 71,95 mm y el mayor máximo absoluto de la muestra, 206,81 mm.Puyo presentó una media de 68,49 mm y un máximo de 159,95 mm, mientras Quito mostró menordispersión, con media de 45,15 mm y máximo de 76,52 mm. La Fig. 2 presenta el caso de Puyo,donde la señal temporal fue más clara; los resultados correspondientes a Guayaquil y Quito se reportancuantitativamente en el párrafo siguiente.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...54
Página 55. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)Fig. 2. Evolución temporal de los máximos anuales de precipitación diaria y tendencia robustade Sen en Puyo.La prueba de Mann-Kendall mostró una señal positiva fuerte en Puyo (τ = 0,281; p < 0,001), conuna pendiente de Sen de 0,405 mm/año e intervalo de confianza positivo. En Guayaquil la pendiente fuetambién positiva, 0,291 mm/año, pero la evidencia quedó próxima al umbral convencional (p = 0,057).Quito no mostró tendencia monotónica (p = 0,869). Estas diferencias confirman que la evaluaciónregional no debe reducirse a una única tendencia nacional y justifican la comparación posterior entreestructuras probabilísticas.B. Selección y parámetros de los modelos GEVLa Tabla 2 resume los AICc de M0-M3 y los BIC relevantes para la decisión final. En Guayaquil,M3 presentó una ventaja amplia respecto de los demás modelos. En Puyo, M1 obtuvo el menor AICc yM3 quedó dentro de dos unidades, pero la incorporación de RONI no produjo una mejora significativaadicional. En Quito, M2 alcanzó el AICc mínimo; sin embargo, M0 quedó dentro de ∆AICc = 1,73 ypresentó un BIC menor que M2 (559,255 frente a 559,604). Dado que ∆AICc = 1,73 < 2, M0 y M2presentan soporte comparable según AICc. En ausencia de una mejora concluyente del modelo completoy al presentar M0 un BIC ligeramente menor, se priorizó la especificación con menos parámetros paraevitar interpretar como no estacionariedad una ganancia marginal de ajuste. Por parsimonia se retuvoM0 como formulación principal y M2 se mantuvo únicamente como análisis de sensibilidad.Tabla 2. Comparación de modelos GEV estacionarios y no estacionarios mediante AICc y BIC.CiudadAICcM0AICcM1AICcM2AICcM3BIC M0 / candidato AICcmínimoModeloretenidoGuayaquil723,715717,660713,806702,727730,330 / 713,444 (M3)M3Puyo685,497680,954686,224682,866692,111 / 689,652 (M1)M1Quito552,640554,577550,906552,498559,255 / 559,604 (M2)M0Las pruebas de razón de verosimilitudes respaldaron M3 frente a M0 en Guayaquil (LR = 25,52;p < 0,001) y M1 frente a M0 en Puyo (LR = 6,78; p = 0,009). En Quito, M2 frente a M0 alcanzó unasignificancia marginal (p = 0,046), pero el modelo completo M3 no mejoró de forma concluyente a M0y la penalización del BIC favoreció la solución estacionaria. Esta decisión conservadora es coherente conla recomendación de no introducir no estacionariedad cuando la ganancia depende de una estructurapoco robusta [5].Los parámetros e intervalos bootstrap de los modelos retenidos se presentan en la Tabla 3. EnGuayaquil, los coeficientes de tiempo y RONI fueron positivos y sus intervalos del 95 % excluyeronPalacios M. et al.Modelización no estacionaria de precipitaciones extremas...55
Página 56. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)cero. En Puyo, el coeficiente temporal también fue positivo con intervalo completamente mayor quecero. Como tiempo y RONI se estandarizaron, estos coeficientes representan cambios asociados a unadesviación estándar de la covariable y no cambios por año o por una unidad RONI. Para Quito, alretener M0, la localización se mantuvo constante. Los diagnósticos bootstrap de Kolmogorov-Smirnov,Cramér-von Mises y Anderson-Darling presentaron valores p superiores a 0,38 en los tres casos, sinevidencia de falta de ajuste.Tabla 3. Parámetros estimados de los modelos GEV seleccionados e intervalos de confianzabootstrap del 95 %.CiudadModeloβ0βtiempoβRONIσξGuayaquilM358,48[53,79; 63,89]8,58[4,51; 12,91]12,30[8,01; 17,06]19,64[15,15; 22,93]0,10[−0,09; 0,32]PuyoM156,19[52,07; 60,83]5,08[1,42; 8,80]—16,87[13,13; 20,07]0,14[−0,05; 0,33]QuitoM040,90[39,11; 42,95]——8,15[6,85; 9,47]−0,06[−0,21; 0,10]La Fig. 3 muestra la evolución de µ en Guayaquil, donde el modelo M3 combina el desplazamientotemporal con las oscilaciones asociadas a RONI. Puyo presentó una evolución aproximadamente li-neal bajo M1, mientras que Quito mantuvo una localización constante bajo M0, comportamientos yacuantificados en la Tabla 3.Fig. 3. Evolución del parámetro de localización del modelo GEV no estacionario seleccionadopara Guayaquil.C. Niveles de retorno y control de robustezLa Tabla 4 presenta los niveles de retorno condicionados. En Guayaquil, manteniendo RONI = 0,el cuantil de 50 años aumentó de 135,9 mm en 1951 a 165,0 mm en 2025, equivalente a una diferenciaaproximada del 21 %. Para 2025, el mismo cuantil varió entre 150,8 mm con RONI = −1 y 179,2mm con RONI = +1. En Puyo, el nivel de 50 años aumentó de 134,4 a 151,6 mm entre 1951 y2025. Quito, bajo M0, mantuvo niveles estacionarios de 58,0; 63,0 y 69,2 mm para 10, 20 y 50 años,respectivamente.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...56
Página 57. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)Tabla 4. Niveles de retorno estimados bajo los escenarios seleccionados e intervalos de confianzadel 95 %.CiudadEscenario10 años (mm)20 años (mm)50 años (mm)Guayaquil1951; RONI = 091,8[76,3; 107,8]110,1[89,5; 135,5]135,9[104,9; 184,0]Guayaquil2025; RONI = −1106,7[90,7; 123,7]125,0[103,8; 151,6]150,8[118,7; 200,7]Guayaquil2025; RONI = 0120,9[105,6; 137,2]139,2[118,5; 165,3]165,0[133,4; 214,6]Guayaquil2025; RONI =+1135,1[118,8; 152,3]153,4[131,8; 179,3]179,2[146,8; 227,6]Puyo195192,0[78,4; 107,4]109,3[90,4; 133,6]134,4[104,3; 178,0]Puyo1988100,6[88,6; 114,9]117,9[99,7; 141,4]143,0[113,0; 185,6]Puyo2025109,2[95,5; 124,4]126,5[107,1; 150,5]151,6[120,4; 196,1]QuitoEstacionario58,0[55,4; 62,1]63,0[59,6; 68,9]69,2[64,0; 78,7]La Fig. 4 concentra la comparación gráfica en Guayaquil, donde la influencia conjunta del tiempo yRONI produjo las diferencias más claras entre niveles de retorno. Los resultados completos de Puyo yQuito se conservan en la Tabla 4. En Guayaquil, el desplazamiento asociado a RONI es comparable enmagnitud al cambio temporal acumulado, lo que concuerda con la fuerte modulación de la precipitacióncostera por ENSO descrita para Ecuador [6], [8]. Los intervalos aumentaron en amplitud conforme seextendió el periodo de retorno, reflejando la mayor extrapolación hacia la cola de la distribución y laincertidumbre de los parámetros GEV. En aplicaciones hidrológicas, los cuantiles de 50 años debeninterpretarse, por tanto, junto con sus intervalos de confianza y no como umbrales deterministas dediseño.Fig. 4. Niveles de retorno estimados para Guayaquil en 2025 bajo RONI = −1, 0 y +1. Losintervalos de confianza del 95 % se presentan en la Tabla 4.El control de robustez con CHIRPS introdujo un matiz relevante. En Guayaquil, ERA5-Land secorrelacionó significativamente con DAILY_RNL (r = 0,496; p = 0,007) y DAILY_SAT (r = 0,509;p = 0,006); las tres series mostraron pendientes positivas en 1998–2025, aunque no significativas en eseintervalo reducido. Por tanto, la dirección general del resultado costero es razonablemente consistentePalacios M. et al.Modelización no estacionaria de precipitaciones extremas...57
Página 58. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)entre productos, aun cuando CHIRPS estimó máximos medios superiores.En Puyo, en cambio, ERA5-Land no mostró correlación con DAILY_RNL ni DAILY_SAT y CHIRPSprodujo máximos considerablemente mayores. La tendencia de ERA5-Land fue significativa en 1998–2025, pero no fue confirmada por ninguna variante CHIRPS. En consecuencia, la selección M1 describede manera válida la serie ERA5-Land de 1951–2025, pero su interpretación climática debe considerarsedependiente del producto. La discrepancia es consistente con diferencias en la construcción de los pro-ductos y en la asignación temporal de la precipitación diaria. ERA5-Land deriva la serie de un reanálisishorario, mientras que las variantes CHIRPS distribuyen los mismos totales pentadales utilizando ERA5en DAILY_RNL e IMERG en DAILY_SAT. En un entorno amazónico con precipitación convectiva in-tensa y elevada variabilidad espacial, estas diferencias pueden modificar tanto la magnitud como el díaasignado a los máximos. Por ello, la señal observada en Puyo se considera sensible al producto y no unatendencia regional confirmada. La discrepancia impide extrapolar este resultado a una intensificaciónamazónica general y refuerza la cautela frente a conclusiones no estacionarias basadas en una solafuente [5], [18].Quito presentó una situación distinta: DAILY_SAT reprodujo prácticamente la media de los máxi-mos de ERA5-Land, aunque la correspondencia interanual fue débil, y ninguna de las tres series mostróuna tendencia significativa robusta en el periodo común. Este resultado respalda la decisión de retenerM0. En conjunto, el control de robustez confirma que la no estacionariedad no puede asumirse comouna propiedad uniforme del territorio y que el producto y la convención temporal de precipitaciónconstituyen fuentes adicionales de incertidumbre al analizar máximos diarios en malla [18].D. Alcances y limitacionesEl alcance de los resultados es puntual y probabilístico: cada localidad se representó mediante unpíxel de productos en malla y no mediante una red de estaciones, por lo que no se infieren cambiospara regiones completas ni se atribuyen causalmente a cambio climático.Además, se permitió noestacionariedad únicamente en el parámetro de localización y el método de máximos anuales reduce lainformación de eventos subanuales. El contraste CHIRPS tampoco constituye validación independiente aescala diaria, porque sus variantes comparten totales pentadales y difieren en la desagregación temporal.Estas limitaciones se incorporaron en la interpretación, especialmente para Puyo.CONCLUSIONESLa comparación de modelos GEV mostró que la estructura probabilística de los máximos diarios deprecipitación difiere entre las tres localidades analizadas. Guayaquil presentó la evidencia más sólidade no estacionariedad, con un parámetro de localización asociado simultáneamente al tiempo y aRONI. Este resultado indica que la variabilidad interanual vinculada a ENSO y el componente temporalaportan información complementaria para caracterizar los extremos costeros, sin implicar por sí solosuna atribución causal a cambio climático.Puyo presentó una señal temporal estadísticamente consistente dentro de ERA5-Land, pero el con-trol con CHIRPS no la corroboró. Por tanto, el resultado debe interpretarse como dependiente delconjunto de datos y no como evidencia universal de intensificación amazónica. Quito permaneció com-patible con estacionariedad cuando se consideraron conjuntamente AICc, BIC, parsimonia y robustez,a pesar de que una formulación con RONI obtuvo una mejora marginal en AICc.Los niveles de retorno condicionados mostraron que asumir estacionariedad puede modificar demanera sustantiva la estimación de eventos raros cuando existe una covariable climática relevante.Sin embargo, los intervalos bootstrap se ampliaron para retornos largos, por lo que la incertidumbredebe acompañar cualquier uso aplicado.La principal implicación metodológica es que la decisiónentre modelos estacionarios y no estacionarios debe combinar evidencia estadística, coherencia física,parsimonia y sensibilidad al producto de precipitación.Como extensión, resulta conveniente incorporar estaciones pluviométricas de largo registro o pro-ductos regionales validados para contrastar los resultados en Puyo y otras zonas amazónicas, así comoevaluar estructuras espaciales y formulaciones alternativas para escala y forma. El flujo de análisisconserva la base maestra, los archivos fuente de tablas y figuras y la configuración del remuestreo, loque facilita su reproducción y adaptación a otras localidades con datos de acceso abierto.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...58
Página 59. Modelización no estacionaria de precipitaciones extremas mediante teoría de valores extremos y datos meteorológicosISSN-e: 2737-6419Período: octubre-diciembre de 2026Revista AtheneaVol.7, Número 26, (pp. 50–59)REFERENCIAS[1] R. W. Katz, M. B. Parlange, and P. Naveau, “Statistics of extremes in hydrology,” Adv. WaterResour., vol. 25, no. 8–12, pp. 1287–1304, Aug. 2002, doi: 10.1016/S0309-1708(02)00056-8.[2] S. Coles, “An introduction to statistical modeling of extreme values,” 2001, doi:10.1007/978-1-4471-3675-0.[3] P. C. D. Milly et al., “Climate change: Stationarity is dead: Whither water management?” Science(1979)., vol. 319, no. 5863, pp. 573–574, Feb. 2008, doi: 10.1126/SCIENCE.1151915.[4] G. Villarini, J. A. Smith, F. Serinaldi, J. Bales, P. D. Bates, and W. F. Krajewski, “Flood frequencyanalysis for nonstationary annual peak records in an urban drainage basin,” Adv. Water Resour.,vol. 32, no. 8, pp. 1255–1266, Aug. 2009, doi: 10.1016/J.ADVWATRES.2009.05.003.[5] F. Serinaldi and C. G. Kilsby, “Stationarity is undead: Uncertainty dominates the distribution ofextremes,” Adv. Water Resour., vol. 77, pp. 17–36, Mar. 2015, doi: 10.1016/J.ADVWATRES.2014.12.013.[6] D. R. Thielen et al., “Effect of extreme El Niño events on the precipitation of Ecuador,” NaturalHazards and Earth System Sciences, vol. 23, no. 4, pp. 1507–1527, Apr. 2023, doi: 10.5194/NHESS-23-1507-2023.[7] G. Urgilés, R. Célleri, J. Bendix, and J. Orellana-Alvear, “Identification of spatio-temporal patternsin extreme rainfall events in the Tropical Andes: A clustering analysis approach,” MeteorologicalApplications, vol. 31, no. 5, p. e70005, Sep. 2024, doi: 10.1002/MET.70005.[8] C. P. Acosta, B. R. de Fonseca, and T. L. Doval, “Modulation of rainfall variability and potentialpredictability in Ecuador: ENSO’s evolving influence (1981–2018),” MAUSAM, vol. 77, no. 4, pp.1313–1346, Oct. 2026, doi: 10.54302/MAUSAM.V77I4.7456.[9] K. P. Burnham and D. R. Anderson, “Model selection and multimodel inference,” 2004, doi:10.1007/B97636.[10] D. Panagoulia, P. Economou, and C. Caroni, “Stationary and nonstationary generalized extremevalue modelling of extreme precipitation over a mountainous area under climate change,” Envi-ronmetrics, vol. 25, no. 1, pp. 29–43, Feb. 2014, doi: 10.1002/ENV.2252.[11] B. Efron and R. J. Tibshirani, “An introduction to the bootstrap,” May 1994, doi: 10.1201/9780429246593.[12] J. D. Salas and J. Obeysekera, “Revisiting the concepts of return period and risk for nonstationaryhydrologic extreme events,” J. Hydrol. Eng., vol. 19, no. 3, pp. 554–568, Apr. 2013, doi: 10.1061/(ASCE)HE.1943-5584.0000820.[13] J. Muñoz-Sabater et al., “ERA5-Land: a state-of-the-art global reanalysis dataset for land ap-plications,” Earth Syst. Sci. Data, vol. 13, no. 9, pp. 4349–4383, Sep. 2021, doi: 10.5194/ESSD-13-4349-2021.[14] N. Gorelick, M. Hancher, M. Dixon, S. Ilyushchenko, D. Thau, and R. Moore, “Google EarthEngine: Planetary-scale geospatial analysis for everyone,” Remote Sens. Environ., vol. 202, pp.18–27, Dec. 2017, doi: 10.1016/J.RSE.2017.06.031.[15] NOAA Climate Prediction Center, “Relative oceanic Niño index (RONI),” accessed: Oct. 21, 2026.[Online]. Available: https://www.cpc.ncep.noaa.gov/products/analysis_monitoring/enso/roni/[16] H. B. Mann, “Nonparametric tests against trend,” Econometrica, vol. 13, no. 3, p. 245, Jul. 1945,doi: 10.2307/1907187.[17] P. K. Sen, “Estimates of the regression coefficient based on Kendall’s tau,” J. Am. Stat. Assoc.,vol. 63, no. 324, pp. 1379–1389, 1968, doi: 10.1080/01621459.1968.10480934.[18] C. Funk et al., “The climate hazards center infrared precipitation with stations, version 3,” Scien-tific Data, vol. 13, no. 1, p. 718, Apr. 2026, doi: 10.1038/s41597-026-07096-4.Palacios M. et al.Modelización no estacionaria de precipitaciones extremas...59