Dos registros pueden contener las mismas cuatro amplitudes y describir evoluciones distintas. En uno las mediciones están separadas por segundos; en el otro, entre la segunda y la tercera pasan varias horas. Numerar las filas como si cada paso durara lo mismo borra parte de esa información, ya que una persistencia elevada a cinco minutos no expresa lo mismo que esa misma persistencia a cinco horas.

Un modelo en tiempo continuo especifica cómo evoluciona el proceso entre los instantes de observación; cada intervalo real determina la transición que corresponde usar. Así se puede trabajar directamente con marcas temporales irregulares, sin completar antes una rejilla por interpolación. Aquí los tiempos se considerarán conocidos y fijados desde fuera del proceso. Si el estado influye en cuándo se mide, el mecanismo de observación requiere un tratamiento aparte.

Que el tiempo sea continuo no obliga a que la variable de estado tome valores continuos. Una cadena puede saltar entre un conjunto finito de estados en cualquier instante, mientras que una difusión puede recorrer de forma continua un espacio de amplitudes. Las cadenas se describen con intensidades de salto y las difusiones con una ecuación diferencial estocástica, pero ambas necesitan traducir sus parámetros a la duración efectiva entre mediciones.

De una tasa de salto a una probabilidad de transición

En una cadena de Markov de tiempo continuo, la matriz generadora QQ contiene intensidades fuera de la diagonal. Para estados distintos ii y jj,

qij=lim⁡δ↓0P(S(t+δ)=j∣S(t)=i)δ.q_{ij}=\lim_{\delta\downarrow0} \frac{P(S(t+\delta)=j\mid S(t)=i)}{\delta}.

Estas cantidades son no negativas y se miden en unidades de tiempo inverso. La diagonal se define por qii=−∑j≠iqijq_{ii}=-\sum_{j\ne i}q_{ij}, con lo que cada fila suma cero. Un generador no es una matriz de probabilidades, ya que sus elementos diagonales pueden ser negativos y los de fuera de la diagonal pueden superar uno si se expresan en una unidad temporal suficientemente grande.

Consideremos dos estados instrumentales, AA y BB, con generador

Q=(−0,40,40,1−0,1),Q= \begin{pmatrix} -0{,}4&0{,}4\\ 0{,}1&-0{,}1 \end{pmatrix},

expresado por hora. El estado AA tiene una intensidad de salida de 0,40{,}4 por hora y BB una de 0,10{,}1 (son valores didácticos). En un intervalo muy breve, la probabilidad de pasar de AA a BB se aproxima por 0,4δ0{,}4\delta. Esa aproximación local no autoriza a multiplicar la tasa por cualquier duración, ya que con cuatro horas daría 1,61{,}6, una probabilidad imposible.

La transición completa de una cadena homogénea se obtiene con la exponencial matricial,

P(Δ)=exp⁡(QΔ)=I+QΔ+(QΔ)22!+⋯ .P(\Delta)=\exp(Q\Delta) =I+Q\Delta+\frac{(Q\Delta)^2}{2!}+\cdots.

La exponencial se aplica a la matriz, a través de sus productos sucesivos. Aplicar la función exponencial elemento por elemento da otra matriz, que en general ni siquiera conserva la suma unitaria de las filas. Para un generador válido de dimensión finita, la exponencial matricial entrega probabilidades no negativas con filas que suman uno.

En nuestro ejemplo la expresión es sencilla,

P(Δ)=(0,2+0,8e−0,5Δ0,8(1−e−0,5Δ)0,2(1−e−0,5Δ)0,8+0,2e−0,5Δ).P(\Delta)= \begin{pmatrix} 0{,}2+0{,}8e^{-0{,}5\Delta}&0{,}8(1-e^{-0{,}5\Delta})\\ 0{,}2(1-e^{-0{,}5\Delta})&0{,}8+0{,}2e^{-0{,}5\Delta} \end{pmatrix}.

Partiendo de AA, la probabilidad de encontrar el sistema en BB es aproximadamente 0,0940 al cabo de un cuarto de hora, 0,3148 al cabo de una hora y 0,6917 al cabo de cuatro. La misma intensidad produce probabilidades distintas porque el tiempo disponible para cambiar es distinto. Para intervalos muy largos, la distribución de destino tiende a 0,2 y 0,8, las ocupaciones estacionarias de esta cadena.

Si el estado se observa exactamente en los tiempos t0,…,tnt_0,\ldots,t_n, la verosimilitud condicionada al estado inicial es el producto

∏k=1nPsk−1,sk(tk−tk−1).\prod_{k=1}^{n} P_{s_{k-1},s_k}(t_k-t_{k-1}).

Cada factor usa su propio intervalo. Si los estados se observan con error o solo a través de una señal, se suma sobre las trayectorias de estados posibles en los instantes de medición, incorporando las distribuciones de emisión. Es la recurrencia de un HMM con una matriz de transición propia para cada intervalo, como la que presenta Jackson (2011) para observaciones de panel junto con su extensión a estados registrados mediante marcadores imperfectos.

La exponencial matricial integra todas las historias posibles entre dos observaciones, también las que contienen varios saltos, por lo que para evaluar la probabilidad de los extremos no hace falta inventar estados intermedios. Sí harían falta cálculos adicionales si se quisiera inferir la trayectoria interior, el número de saltos o el tiempo acumulado en cada estado.

La construcción llamada uniformización hace explícita esa suma. Se elige una tasa común al menos tan grande como todas las intensidades de salida y se generan oportunidades de transición mediante un proceso de Poisson; en cada oportunidad se aplica una matriz discreta que admite también conservar el estado. En nuestro ejemplo basta la tasa 0,40{,}4, con matriz discreta R=I+Q/0,4R=I+Q/0{,}4, cuyas filas son (0; 1)(0;\,1) y (0,25; 0,75)(0{,}25;\,0{,}75). El número de oportunidades durante cuatro horas tiene media 1,61{,}6.

Dado ese número, se eleva RR a la potencia correspondiente y después se ponderan las potencias por las probabilidades de Poisson, con lo que se recupera exactamente la exponencial matricial. Las oportunidades que conservan el estado son eventos virtuales del procedimiento y no corresponden a saltos físicos adicionales. Gracias a esa separación pueden calcularse probabilidades y simularse trayectorias sin confundir el reloj auxiliar con los cambios reales del sistema. Hobolth y Jensen (2011) usan esta construcción para estudiar estadísticas condicionadas a los extremos.

Permanencia, retorno y calendario de visitas

Al entrar en AA, el tiempo hasta el primer salto sigue una distribución exponencial de tasa 0,4, con una media de 2,5 horas y función de supervivencia

P(Tsalida>u∣S(0)=A)=e−0,4u.P(T_{\mathrm{salida}}>u\mid S(0)=A)=e^{-0{,}4u}.

Para BB, la media hasta el primer salto es de diez horas. Estas duraciones físicas suponen una trayectoria observada con el detalle suficiente para reconocer su primer cambio de estado.

Al cabo de cuatro horas, la probabilidad de no haber salido nunca de AA es aproximadamente 0,20190{,}2019, mientras que la de encontrar el sistema en AA al final es mayor, PAA(4)≈0,3083P_{AA}(4)\approx0{,}3083. La diferencia, 0,10640{,}1064, corresponde a trayectorias que salieron y volvieron. Si se sabe que ambos extremos son AA, la probabilidad condicional de que haya habido al menos una salida y un regreso es cercana a 0,34510{,}3451.

Encontrar al final el mismo estado del principio deja, por tanto, incertidumbre sobre la historia intermedia. Hobolth y Jensen (2011) desarrollan el cálculo de estadísticas de trayectoria condicionadas a sus extremos, como tiempos de ocupación y números de saltos. Esas expectativas combinan la información de los extremos con la del generador, pero no hacen observables los eventos intermedios ni determinan sus instantes exactos.

La diferencia reaparece al contar rachas en una rejilla regular. Si se parte observando AA y se repite la observación cada hora, la probabilidad de encontrar BB en la visita siguiente es 0,3148. El número de intervalos hasta la primera visita que registra BB tiene una distribución geométrica, con un tiempo medio de aproximadamente 3,1769 horas, cuando el tiempo medio físico hasta el primer salto sigue siendo de 2,5. La visita puede detectar el cambio con retraso o pasar por alto una excursión completa a BB.

Con visitas cada dos horas, el tiempo medio hasta el primer registro de BB sube a aproximadamente 3,9549 horas. Esa cantidad depende de la rejilla de medición y no refleja un cambio de la tasa física del proceso, que se ha mantenido constante. La media se obtiene dividiendo dos horas por PAB(2)P_{AB}(2), sin convertir directamente la duración media de la cadena continua en un número entero de visitas.

Con intervalos irregulares ni siquiera existe una única probabilidad de éxito para cada nueva visita. Si las observaciones anteriores han seguido registrando AA, la probabilidad de volver a encontrarlo es PAA(Δk)P_{AA}(\Delta_k) y la de una racha concreta depende del producto de esos factores. Ajustar una geométrica al número de lecturas consecutivas con la misma etiqueta mezcla entonces duraciones físicas con el calendario de observación.

La distinción ayuda a decidir qué resultado comunicar. En un proceso de alternancia quizá interese el tiempo hasta salir de un estado; para organizar inspecciones, el tiempo hasta detectar otro estado; para cuantificar una exposición, el total de horas acumuladas en AA, sumando varias visitas a ese estado dentro del mismo intervalo. El generador define las tres cantidades, aunque los cálculos y la información observada que las sostienen son diferentes.

Composición de intervalos y orden de los regímenes

La transición homogénea cumple P(a+b)=P(a)P(b).P(a+b)=P(a)P(b).

El producto suma sobre los estados intermedios posibles. Si no hubo ninguna medición entre la primera hora y las tres siguientes, propagar la distribución en esos dos tramos da lo mismo que propagarla directamente durante cuatro horas. Partir un intervalo para facilitar el cálculo no añade información.

Cuando sí hay una medición intermedia la situación cambia. Se propaga hasta ese instante, se incorpora su emisión y se vuelve a propagar. El término de emisión se interpone entre las dos matrices y modifica los pesos de las trayectorias. La igualdad anterior vale para la propagación sin información nueva y no permite eliminar una medición real situada en medio del intervalo.

Si el generador cambia con el tiempo, importa además el orden de los tramos. Supongamos que durante una hora rige QQ y durante la siguiente rige

Q2=(−0,10,10,4−0,4).Q_2= \begin{pmatrix} -0{,}1&0{,}1\\ 0{,}4&-0{,}4 \end{pmatrix}.

Con las distribuciones escritas como vectores fila, la transición total es exp⁡(Q)exp⁡(Q2)\exp(Q)\exp(Q_2). Partiendo de AA, la probabilidad final de BB es aproximadamente 0,2696; si se invierte el orden de los generadores, pasa a 0,3625. Las matrices no conmutan y el régimen que actúa después recibe una distribución distinta según lo que ocurrió antes. Reemplazar los dos generadores por una tasa promedio borra esa secuencia y en general no reproduce la transición.

Cuando los parámetros cambian de forma continua, la propagación se obtiene resolviendo las ecuaciones de evolución correspondientes o con aproximaciones numéricas controladas. Una expresión como la exponencial de la integral del generador solo es válida bajo condiciones adicionales de conmutación. La representación por tramos constantes es transparente cuando se conoce el momento de un cambio externo y cada tramo tiene parámetros propios.

El camino inverso, desde una matriz discreta ajustada hasta un generador, no está asegurado. La matriz

(0,20,80,70,3)\begin{pmatrix} 0{,}2&0{,}8\\ 0{,}7&0{,}3 \end{pmatrix}

es una matriz de transición discreta válida, pero tiene determinante −0,5-0{,}5. No puede ser la exponencial de una matriz real, porque el determinante de esa exponencial es la exponencial de la traza de la matriz y, por tanto, positivo. En particular, no puede representar el muestreo de una cadena homogénea continua sobre esos mismos dos estados. Ajustar un generador impone una familia de transiciones compatibles entre intervalos; calcular un logaritmo matricial sin comprobar sus propiedades no resuelve esta incompatibilidad.

Una difusión con retorno gradual

Dos curvas temporales comparan media y varianza de la transición exacta y un paso Euler para intervalos de cero a tres horas.
Figura 326.1. La difusión Ornstein–Uhlenbeck exacta tiene media 10,8925 y varianza 3,8009 tras tres horas, mientras un solo paso Euler produce media 8 y varianza 12. Ampliar SVG ↗

Para una variable cuya amplitud cambia de forma continua, consideremos el proceso

dXt=−κ(Xt−μ) dt+σ dWt.dX_t=-\kappa(X_t-\mu)\,dt+\sigma\,dW_t.

Aquí μ\mu es un nivel de referencia, κ>0\kappa>0 fija la velocidad de retorno, WtW_t es un movimiento browniano estándar y σ\sigma mide la intensidad de las perturbaciones. Este modelo de Ornstein-Uhlenbeck nació de la descripción del movimiento browniano con relajación (Uhlenbeck y Ornstein, 1930).

Si el tiempo se mide en horas y la amplitud en unidades instrumentales, κ\kappa se expresa por hora y σ\sigma en unidades instrumentales por raíz de hora. El incremento browniano durante un intervalo tiene varianza igual a la longitud del intervalo. Su desviación típica crece con la raíz cuadrada del tiempo, mientras que el término de deriva se acumula al principio en proporción al tiempo.

La notación diferencial resume una ecuación integral y no supone que una trayectoria browniana tenga una derivada ordinaria que pueda evaluarse como una señal suave. En este modelo el ruido es aditivo y las soluciones tienen trayectorias continuas, aunque irregulares a escala local, que una línea trazada entre observaciones no reproduce.

La ecuación lineal puede resolverse exactamente entre dos tiempos,

Xt+Δ=μ+e−κΔ(Xt−μ)+ηΔ,X_{t+\Delta} =\mu+e^{-\kappa\Delta}(X_t-\mu)+\eta_\Delta,

donde la perturbación es normal, independiente del estado inicial y de las perturbaciones de intervalos anteriores que no se solapan, con varianza

V(Δ)=σ22κ(1−e−2κΔ).V(\Delta)= \frac{\sigma^2}{2\kappa} \left(1-e^{-2\kappa\Delta}\right).

La normalidad de esta transición se debe a que se integra una función determinista contra el movimiento browniano. No requiere que la distribución inicial de XtX_t sea normal para que la ley condicional, dado su valor, tenga esta forma. Si se propaga una distribución inicial no gaussiana, su transformación más el ruido añadido pueden dar una distribución marginal también no gaussiana.

El adjetivo «exacta» indica que los estados en los tiempos elegidos pueden generarse con su distribución correcta sin reconstruir los intermedios del movimiento browniano. Si la pregunta depende de lo ocurrido entre dos mediciones (si se superó un umbral y luego se regresó, por ejemplo), los extremos no son suficientes. Una trayectoria interpolada linealmente supondría una ausencia de excursiones que el modelo no impone. La distribución de trayectorias condicionadas a sus extremos puede estudiarse, pero es una inferencia adicional que conserva variabilidad entre las trayectorias compatibles.

Con μ=10\mu=10, κ=0,5\kappa=0{,}5 por hora y σ=2\sigma=2, supongamos conocido exactamente Xt=14X_t=14. La tabla recoge las transiciones exactas.

Tiempo transcurrido, en horasMedia de Xt+ΔX_{t+\Delta}Varianza de la transición
0,2513,53000,8848
112,42612,5285
310,89253,8009
510,32833,9730

Desde catorce, la media predicha se acerca a diez sin cruzar ese nivel, aunque una trayectoria individual puede cruzarlo por efecto del ruido. La varianza condicional parte de cero y tiende a cuatro porque el retorno amortigua las perturbaciones antiguas. Suponer que la incertidumbre crece siempre en proporción al tiempo ignoraría esa disipación.

El coeficiente e−κΔe^{-\kappa\Delta} depende del intervalo y vale cerca de 0,6065 a una hora y 0,2231 a tres. Tratar ambos valores como dos estimaciones del mismo coeficiente autorregresivo discreto atribuiría al sistema una diferencia que en realidad crea el calendario. Driver y colaboradores (2017) usan modelos de tiempo continuo para representar relaciones dinámicas con intervalos variables entre mediciones. Con los parámetros continuos de esos modelos pueden expresarse después los efectos para un horizonte temporal común.

La varianza de la perturbación discreta también debe cambiar con el intervalo, porque mantener el coeficiente correcto pero reutilizar una varianza fija para todos los intervalos produciría distribuciones predictivas incorrectas. Para un intervalo muy corto, la expansión de la expresión exacta da V(Δ)≈σ2ΔV(\Delta)\approx\sigma^2\Delta; para uno largo, el mismo término se satura. La aproximación de intervalo corto y el comportamiento asintótico son dos tramos de una misma función.

La composición conserva ambos elementos. Si a(Δ)=e−κΔa(\Delta)=e^{-\kappa\Delta}, entonces

V(Δ1+Δ2)=a(Δ2)2V(Δ1)+V(Δ2).V(\Delta_1+\Delta_2) =a(\Delta_2)^2V(\Delta_1)+V(\Delta_2).

El ruido del primer tramo atraviesa la dinámica del segundo antes de llegar al final y por eso su varianza se multiplica por el cuadrado del coeficiente de propagación. Sumar sin más las dos varianzas daría al ruido antiguo un peso que ya no conserva. La derivación general de estas transiciones lineales se encuentra en Särkkä y Solin (2019).

Transmisión del ruido entre componentes

En varias dimensiones, una ecuación lineal centrada puede escribirse como

dXt=FXt dt+L dWt.d\mathbf X_t=F\mathbf X_t\,dt+L\,d\mathbf W_t.

La matriz FF relaciona las componentes del estado y LL reparte entre ellas un movimiento browniano vectorial estándar. Con coeficientes constantes, la transición exacta tiene matriz

A(Δ)=eFΔA(\Delta)=e^{F\Delta}

y covarianza de perturbación

Ω(Δ)=∫0ΔeFuLLTeFTu du.\Omega(\Delta)= \int_0^\Delta e^{Fu}LL^\mathsf{T}e^{F^\mathsf{T}u}\,du.

La integral expresa que una perturbación aplicada al comienzo puede desplazarse entre componentes antes de la medición siguiente. Aunque el ruido entre directamente en una sola componente, más tarde puede afectar a varias, por lo que la matriz de covarianza discreta no tiene por qué conservar el patrón de ceros de la matriz instantánea.

Un ejemplo muestra el efecto sin complicar la dinámica. Sea XtX_t una posición y VtV_t una velocidad, con dXt=Vt dtdX_t=V_t\,dt y dVt=dWtdV_t=dW_t, en unidades compatibles y con intensidad de ruido igual a uno. La velocidad recibe las perturbaciones y la posición integra la velocidad. La matriz de transición y la covarianza añadida son

A(Δ)=(1Δ01),Ω(Δ)=(Δ3/3Δ2/2Δ2/2Δ).A(\Delta)= \begin{pmatrix}1&\Delta\\0&1\end{pmatrix},\qquad \Omega(\Delta)= \begin{pmatrix} \Delta^3/3&\Delta^2/2\\ \Delta^2/2&\Delta \end{pmatrix}.

Al cabo de dos unidades de tiempo, la varianza adicional de la posición es 8/38/3, la de la velocidad es dos y la covarianza entre ambas también es dos. Una perturbación que eleva la velocidad tiende a elevar la posición posterior y la covarianza registra esa historia compartida.

Aquí I+FΔI+F\Delta coincide exactamente con la exponencial porque F2=0F^2=0. Aun así, tomar LLTΔLL^\mathsf{T}\Delta como covarianza daría varianza adicional nula para la posición. Acertar con la propagación de la media no basta para acertar con la incertidumbre, que exige integrar también cómo se transmite el ruido durante el intervalo.

Para propagar una covarianza previa CC se calcula A(Δ)CA(Δ)T+Ω(Δ)A(\Delta)CA(\Delta)^\mathsf{T}+\Omega(\Delta), cuyo primer término transporta la incertidumbre que ya existía mientras el segundo incorpora las perturbaciones nuevas. Con observaciones irregulares ambos dependen del tiempo transcurrido, aunque el instrumento mantenga la misma precisión. Esa separación ayuda a reconocer si una discrepancia viene del proceso o del dispositivo.

La forma integral permite verificar propiedades del resultado. En cualquier dirección del espacio de estados la varianza añadida es una integral de cuadrados, por lo que no es negativa. La covarianza puede ser singular si hay direcciones que el ruido nunca alcanza. En el ejemplo de posición y velocidad ambas direcciones adquieren incertidumbre y la matriz es definida positiva para cualquier intervalo positivo, pese a que el ruido instantáneo entra solo por la velocidad. Qué direcciones pueden variar depende a la vez de la duración del intervalo y de la transmisión entre componentes.

La aproximación de Euler-Maruyama

Cuando no se dispone de una transición exacta, una aproximación frecuente es el método de Euler-Maruyama. Aplicado al modelo escalar anterior con un paso de longitud hh, da

Xt+h≈Xt−κ(Xt−μ)h+σh Z,Z∼N(0; 1).X_{t+h}\approx X_t-\kappa(X_t-\mu)h+\sigma\sqrt h\,Z, \qquad Z\sim N(0;\,1).

El método evalúa la deriva al comienzo del paso y la mantiene fija durante todo el intervalo. Higham (2001) presenta esta construcción y distingue entre aproximar trayectorias y aproximar distribuciones o expectativas. Sus propiedades de convergencia exigen condiciones sobre los coeficientes y pasos suficientemente pequeños.

En el ejemplo, un único paso de un cuarto de hora da media 13,5 y varianza uno, cerca de los valores exactos 13,5300 y 0,8848. Con un paso de tres horas la media cae a ocho y la varianza sube a doce, cuando la transición exacta daba 10,8925 y 3,8009. La aproximación cruzó el nivel de retorno y multiplicó la incertidumbre.

El cambio de signo aparece porque Euler reemplaza e−κhe^{-\kappa h} por 1−κh1-\kappa h, que con tres horas vale −0,5-0{,}5, cuando el exacto sigue siendo positivo. Si el modelo aproximado se usara para inferir la dinámica a partir de esos intervalos, podría introducir una alternancia que la transición exacta no contiene.

La estabilidad también se altera. Con pasos iguales repetidos, la parte autorregresiva de Euler exige ∣1−κh∣<1\lvert1-\kappa h\rvert<1 para que la varianza estacionaria sea finita con este ruido aditivo. Para κ=0,5\kappa=0{,}5 el paso debe ser menor de cuatro horas; con cinco horas el coeficiente vale −1,5-1{,}5 y las desviaciones se amplifican. El proceso continuo, en cambio, es estable porque su coeficiente exacto está entre cero y uno para cualquier intervalo positivo.

Incluso dentro del rango estable la variabilidad puede deformarse. Con pasos de una hora, la varianza estacionaria de Euler es de aproximadamente 5,3333, frente a cuatro en el modelo continuo. Con pasos de tres horas sube a dieciséis. Por tanto, comprobar solo que una simulación no explote no garantiza que reproduzca la distribución de interés.

El intervalo entre observaciones y el paso de integración son decisiones distintas. Un hueco de tres horas puede recorrerse con muchos pasos internos aunque no haya ninguna medición en ese lapso. Con diez pasos de Euler, partiendo de catorce, la media final es 10,7875 y la varianza 4,1567; con cien pasos, 10,8824 y 3,8341, ya cerca de los valores exactos. Esos estados internos pertenecen al cálculo y no deben presentarse como observaciones recuperadas.

En una simulación aproximada, reducir el paso muestra cuánto cambian las cantidades que se usarán después. Si el objetivo es la probabilidad de superar un umbral, comprobar solo la media puede dejar sin examinar el error que importa. En modelos no lineales o con ruido dependiente del estado el comportamiento puede ser más difícil que en este ejemplo. Disponer de una solución lineal exacta da un caso de comparación en el que el error numérico puede localizarse sin confundirlo con incertidumbre estadística.

Hay dos tipos de comparación numérica. Para saber si una trayectoria aproximada se acerca a la que corresponde a una realización concreta del ruido, una comprobación fuerte acopla ambas a los mismos incrementos brownianos. En una comprobación débil se comparan distribuciones de resultados para estudiar una media, una varianza u otra expectativa. Dos simulaciones con perturbaciones independientes pueden separarse mucho aunque ambas sean exactas. Su diferencia visible no mide por sí sola el error del integrador, del mismo modo que dos trayectorias parecidas no garantizan una covarianza correcta.

Información contenida en intervalos diversos

El generador de una cadena puede estimarse buscando las intensidades que dan máxima verosimilitud a las transiciones observadas, cada una evaluada en su intervalo. Si un estado se vigila de forma continua, se conocen los saltos y los tiempos de exposición, mientras que con observaciones de panel esas cantidades quedan ocultas entre visitas. Reemplazarlas por el número de cambios visibles y por duraciones asignadas a partir de las etiquetas puede alterar la estimación. El cálculo con probabilidades de transición integra los eventos que el registro no distingue.

En la difusión, los intervalos cortos y los largos informan sobre aspectos diferentes. A intervalos cortos, los incrementos dependen todavía mucho del valor inicial y la varianza nueva empieza creciendo casi linealmente. A intervalos largos, la influencia del valor anterior se debilita y la distribución se acerca al régimen estacionario. Si todas las separaciones son muy largas, las observaciones pueden informar bien sobre el nivel medio y la variabilidad estacionaria mientras aportan poca información sobre la rapidez de retorno.

En el modelo escalar esa variabilidad es el cociente σ2/(2κ)\sigma^2/(2\kappa). Si la intensidad de difusión y la velocidad de retorno cambian juntas de modo que ese cociente se mantenga constante, pueden producirse distribuciones marginales similares. Las transiciones a escalas donde todavía queda memoria ayudan a distinguir esos parámetros. La irregularidad aporta varios horizontes, aunque su utilidad depende de cuáles sean y de cuántas observaciones caigan en cada rango. Añadir una separación enorme a un registro hecho de separaciones enormes puede cambiar poco la información sobre la dinámica rápida.

Cuando hay error de medición, las observaciones cercanas cumplen otra función. Muestran cómo varían las diferencias registradas a separaciones pequeñas, donde el proceso todavía ha cambiado poco pero el dispositivo ha vuelto a medir. Con errores independientes cada lectura añade su propia incertidumbre y con errores correlacionados el patrón es otro. Estimar la difusión atribuyendo todo el movimiento observado al proceso puede sobrevalorar sus perturbaciones. La ecuación de observación debe acompañar a la dinámica también cuando los tiempos se conocen con precisión.

La condición inicial forma parte también de la verosimilitud. Un experimento que fija el estado al comenzar permite condicionar a ese valor. Un registro abierto en un momento cualquiera necesita una distribución inicial, que puede ser la estacionaria si esa hipótesis está justificada. Usar la misma transición continua en ambos casos no hace idénticas sus verosimilitudes. Las primeras observaciones informan sobre la inicialización, sobre todo cuando el retorno es lento y su influencia dura una fracción apreciable del registro.

Incorporar una nueva lectura después de un hueco

Supongamos ahora que el estado inicial se conoce con incertidumbre y que su distribución filtrada es normal, de media catorce y varianza uno. El instrumento mide Y=X+εY=X+\varepsilon, con error normal independiente de varianza uno. Se mantienen los parámetros continuos del ejemplo. Antes de observar de nuevo, la media se propaga igual que antes, pero la varianza predictiva del estado incluye además la incertidumbre inicial,

C−(Δ)=e−Δ+4(1−e−Δ)=4−3e−Δ.C^-(\Delta)=e^{-\Delta}+4(1-e^{-\Delta}) =4-3e^{-\Delta}.

Al cabo de un cuarto de hora vale aproximadamente 1,6636; al cabo de tres horas, 3,8506. En ambos casos la varianza predictiva de la lectura añade uno por el error de medición. Importa hacer esa suma después de propagar el estado, ya que el ruido de un instrumento que todavía no ha medido no debe atravesar ficticiamente la dinámica del intervalo anterior.

Que la incertidumbre crezca con el hueco depende aquí de que la varianza inicial sea menor que la estacionaria. Con una varianza inicial de diez, la predicción daría 4+6e−Δ4+6e^{-\Delta}, que decrece hacia cuatro, porque el retorno contrae parte de la dispersión inicial mientras incorpora ruido nuevo. Por tanto, la longitud del hueco no determina por sí sola si la varianza total aumenta o disminuye. Hay que comparar lo que el modelo propaga desde la distribución anterior con lo que añade durante el intervalo.

Si la nueva lectura es doce, la ganancia de la actualización gaussiana es aproximadamente 0,6246 tras el intervalo corto y 0,7938 tras el largo, con medias actualizadas de 12,5744 y 11,7717. En el primer caso la predicción seguía cerca de catorce y la lectura la rebaja; en el segundo la predicción ya había vuelto hacia diez y la lectura la eleva. La misma amplitud observada tiene efectos distintos porque llega después de evoluciones distintas.

El cálculo usa la actualización lineal gaussiana habitual. Lo que aporta el modelo continuo son la media y la covarianza previas adecuadas a cada hueco, sin necesidad de interpolar una serie regular antes de actualizar. Tampoco hay razón para asignar la misma ganancia a todas las observaciones cuando los intervalos que las preceden son diferentes.

La unidad de tiempo debe mantenerse coherente. Para expresar el mismo proceso en minutos, la intensidad de retorno pasa de 0,50{,}5 por hora a 1/1201/120 por minuto y la amplitud de difusión pasa de dos a 2/602/\sqrt{60}, aproximadamente 0,25820{,}2582, por raíz de minuto. Dividir también la difusión por sesenta describiría otro proceso. En la cadena de saltos todas las intensidades del generador se dividen por sesenta, de modo que una transición de sesenta minutos coincide con la que antes se calculaba para una hora.

Las marcas temporales delimitan además qué mide el dispositivo. Una lectura instantánea se relaciona con XtX_t, mientras que un promedio acumulado durante una ventana se relaciona con una integral del estado. Dar a ese promedio la ecuación de una medición puntual altera su distribución y su relación con las demás lecturas. Si la ventana es grande frente a la escala de cambio, puede ampliarse el estado para incluir la integral o derivarse una ecuación de observación apropiada.

Dos ventanas que se solapan comparten un tramo del proceso integrado y sus promedios pueden estar muy relacionados aunque sus marcas temporales finales sean distintas. Esa relación no se explica necesariamente aumentando la persistencia del estado instantáneo, ya que también interviene cuánto tiempo del registro comparte cada par de ventanas. Guardar el comienzo y el final de cada integración, junto con la hora de entrega del dato, permite formular esa parte del modelo.

Algo parecido ocurre cuando cambia la frecuencia de un dispositivo que calcula sus salidas con promedios internos. Obtener más filas no implica observar el estado con la misma precisión durante intervalos más cortos, porque puede haberse modificado el periodo de integración, el error de la medida o la dependencia entre salidas. El calendario externo y la construcción de cada lectura describen operaciones diferentes. Para comparar dos registros conviene expresar ambos con una transición en tiempo físico y con una ecuación de observación acorde con su procedimiento de adquisición.

Consideremos, para terminar, un registro didáctico sin error de medición que aísla la contribución del calendario a la verosimilitud. Los estados 14; 13,8; 11,914;\ 13{,}8;\ 11{,}9 y 10,610{,}6 se registran en el instante cero, al cuarto de hora, a la hora y cuarto y a las cuatro horas y cuarto. Se condiciona al primer valor. Los intervalos que entran en las transiciones son 0,250{,}25, uno y tres horas.

La primera predicción tiene media 13,530013{,}5300 y varianza 0,88480{,}8848. Para la siguiente se parte del estado ya observado, 13,813{,}8, que propagado durante una hora da media 12,304812{,}3048 y varianza 2,52852{,}5285. La última parte de 11,911{,}9 y, al cabo de tres horas, da media 10,423910{,}4239 y varianza 3,80093{,}8009. Las log-densidades condicionales suman aproximadamente −3,9047-3{,}9047.

Cada observación nueva ha pasado a ser el punto de partida de la transición siguiente, mientras su marca temporal ha determinado cuánto se atenúa ese punto de partida y cuánta incertidumbre se agrega. Si la última medición hubiera llegado antes, habría que recalcular ambos términos para un intervalo más corto. Con las marcas registradas, el último pronóstico resulta de propagar durante tres horas el valor 11,911{,}9, aunque esas dos mediciones ocupen filas consecutivas del archivo.

Referencias

Driver, C. C., Oud, J. H. L., & Voelkle, M. C. (2017). Continuous time structural equation modeling with R package ctsem. Journal of Statistical Software, 77(5), 1–35. https://doi.org/10.18637/jss.v077.i05

Higham, D. J. (2001). An algorithmic introduction to numerical simulation of stochastic differential equations. SIAM Review, 43(3), 525–546. https://doi.org/10.1137/S0036144500378302

Hobolth, A., & Jensen, J. L. (2011). Summary statistics for endpoint-conditioned continuous-time Markov chains. Journal of Applied Probability, 48(4), 911–924. https://doi.org/10.1239/jap/1324046009

Jackson, C. H. (2011). Multi-state models for panel data: The msm package for R. Journal of Statistical Software, 38(8), 1–28. https://doi.org/10.18637/jss.v038.i08

Särkkä, S., & Solin, A. (2019). Applied stochastic differential equations. Cambridge University Press. https://doi.org/10.1017/9781108186735

Uhlenbeck, G. E., & Ornstein, L. S. (1930). On the theory of the Brownian motion. Physical Review, 36(5), 823–841. https://doi.org/10.1103/PhysRev.36.823