Correlación con Python – Ejercicio Estación Meteorológica

Referencia: Ross 8.2 p334

De una estación meteorológica se obtiene un archivo.csv con los datos de los sensores disponibles.

2021_10_EstMeteorologica.csv

 

Nota: las instrucciones para la lectura del archivo se describen en:

Archivos.csv con Python – Ejercicio con gráfica de temperatura y Humedad

Para observar la correlación entre las medidas obtenidas de varios sensores, se considera usar las variables:

sensor_EM = ["TEMP","Humidity","Solar_rad","Bar_press."]

La correlación entre las variables se realiza seleccionando las columnas desde una lista, para luego realizar la operación desde una de las funciones de Pandas-Python

# Matriz de Correlación
mediciones = tablaEM[sensor_EM].copy()
correlacion_matriz = mediciones.corr()

obteniendo como resultados:

Matriz de correlación:
                TEMP  Humidity  Solar_rad  Bar_press.
TEMP        1.000000 -0.957730   0.743379   -0.387187
Humidity   -0.957730  1.000000  -0.765766    0.447845
Solar_rad   0.743379 -0.765766   1.000000   -0.093236
Bar_press. -0.387187  0.447845  -0.093236    1.000000

Del resultado se observa que el coeficiente de correlación entre las mediciones de «TEMP» y «Humidity» es alta. Para observar el caso se usa la gráfica entre las dos variables comparando los datos en el tiempo por fechas.
Cuando la temperatura aumenta la Humedad disminuye, lo opuesto tambien es observable.

Comparando solo entre las dos variables:

Para el caso de radiación solar y presión baromérica, el coeficiente de correlación es «bajo».

comparando entre las dos variables solamente:

Instrucciones en Python

# Archivos de estación meteorologica
# correlacion entre sensor_EM
import numpy as np
import pandas as pd
import datetime as dt
import os
import matplotlib.pyplot as plt
from matplotlib.dates import DateFormatter, DayLocator

# INGRESO
archivoEM = "2021_10_EstMeteorologica.csv"
sensor_EM = ["TEMP","Humidity","Solar_rad","Bar_press."]
sensor_EM_u = ["(°C)","(%)","(kW)","(hPa)"]

# PROCEDIMIENTO
# tabla estacion meteorologica
tablaEM = pd.read_csv(archivoEM, sep=';',decimal=',')
tablaEM = pd.DataFrame(tablaEM)
tablaEM_n = len(tablaEM)

# Matriz de Correlación
mediciones = tablaEM[sensor_EM].copy()
correlacion_matriz = mediciones.corr()

# SALIDA
print('Matriz de correlación:')
print(correlacion_matriz)

Para añadir las gráficas, se añaden instrucciones para seleccionar las variables en intervalos de fechas.

# tabla estacion meteorologica
fechaformatoEM = "%d/%m/%Y %H:%M:%S"
# fechas concatenando columnas de texto
tablaEM['fecha'] = tablaEM['Date']+' '+tablaEM['Time']
tablaEM['fecha'] = pd.to_datetime(tablaEM['fecha'],
                                  format=fechaformatoEM)
# grafica
sensorA = "TEMP"
sensorB = "Humidity"
fechaA = "2021/10/01"
fechaB = "2021/10/07"

ti = tablaEM["fecha"].loc[tablaEM["fecha"].between(fechaA,fechaB)]
xi = tablaEM[sensorA].loc[tablaEM["fecha"].between(fechaA,fechaB)]
yi = tablaEM[sensorB].loc[tablaEM["fecha"].between(fechaA,fechaB)]

# grafica compara sensores en tiempo
fig_Pt, graf1 = plt.subplots()
graf1.scatter(ti,xi, marker = '.',
              label=sensorA)

eje12 = graf1.twinx()  # segundo eje del grafico
eje12.scatter(ti,yi, marker = '.',
              color='orange',label=sensorB)
graf1.set_ylabel(sensorA, color='blue')
eje12.set_ylabel(sensorB, color='orange')
plt.show()

# grafica compara sensores
plt.scatter(xi,yi)
plt.xlabel(sensorA, color='blue')
plt.ylabel(sensorB, color='orange')
plt.show()

 

Correlación bits con ruido

Dada una secuencia de bits, al transmitirlos se distorsionó con ruido aditivo. En el receptor será dificil discriminar los valores de los bits tomados de la señal recibida.

La correlación permite «limpiar» un poco la señal en el receptor y una mejor estimación de los bits recibidos.

Para el ejemplo se usa una secuencia de bits 01101001 como señal a transmitir, el ruido aditivo es de tipo normal con media m=0 y varianza σ2=1

Para referencia del proceso, se usa un punto rojo en el centro de la ventana de tiempo de cada bit.

La correlación se realiza con un bit y la señal con ruido.

Observe que luego de la correlación es más sencillo discriminar si lo recibido fué un bit 0 o un bit 1.

señal de bits sin ruido
[0 1 1 0 1 0 0 1]
estimado de bits en el receptor
[ 0.002  0.985  1.07   0.054  0.947 -0.079 -0.006  1.149]

Instrucciones en Python

# secuencia de bits con Ruido
# usar correlación para mejorar la senal
# estimar el bit recibido
# Tarea: convertir estimado en bits

import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as st

# INGRESO

senalbit = np.array([0, 1, 1, 0, 1, 0, 0, 1])

# PROCEDIMIENTO

# señal moduladora son escalones por cada bit
# o forma rectangular de duración anchobit
anchobit = 128
mitadbit = anchobit//2
senal = np.repeat(senalbit, anchobit)
muestras=len(senal)

# reloj para observar la mitad del bit
reloj = np.arange(mitadbit, muestras, anchobit)

# Ruido normal o Gausiano
media = 0
varianza = 1
ruido = st.norm.rvs(media,varianza,muestras)

# Añade el ruido a la señal
senalruido = senal + ruido

# Referencia de un bit para correlación
modelo = np.ones(anchobit)

# correlación modelo de bits
correlacion = np.correlate(senalruido, modelo, mode='same')
# normaliza el resultado [0,1]
correlacion = correlacion / anchobit

# estimado de bits en el receptor
estimado=correlacion[reloj]

# SALIDA
print('señal de bits sin ruido')
print(senalbit)
print('estimado de bits en el receptor')
np.set_printoptions(precision=3)
print(estimado)

#GRAFICA
plt.figure(1)       # define la grafica
plt.suptitle('señal binaria con ruido gaussiano')

plt.subplot(311)
plt.plot(senal)
plt.plot(reloj, senal[reloj], 'ro')
plt.ylabel('señal original')
plt.margins(0,0.05)

plt.subplot(312)
plt.plot(senalruido)
plt.plot(reloj, senalruido[reloj], 'ro')
plt.ylabel('senal con ruido')
plt.margins(0,0.05)

plt.subplot(313)
plt.plot(correlacion)
plt.ylabel('correlación')
plt.plot(reloj, correlacion[reloj], 'ro')
plt.margins(0,0.05)

plt.show()

Autocorrelación Propiedades

Referencia León-García 9.6 p522

Potencia Promedio

La función autocorrelación a τ=0 entrega la potencia promedio (segundo momento) del proceso:

R_{X}(0) = E[X(t)^2]

para todo valor de t

Función par respecto a τ

R_{X}(\tau) = E[X(t+\tau)X(t)] = E[X(t)X(t+\tau)] = R_{X}(-\tau)

Mide la tasa de cambio

La función de autocorrelación es una medida de la tasa de cambio de un proceso aleatorio

P[|X(t+\tau) - X(t)| > \epsilon] = P[(X(t+\tau) - X(t))^2 > \epsilon ^2] \leq \frac{E[(X(t+\tau) - X(t))^2]}{\epsilon ^2} =\frac{2\{R_X(0) - R_X(\tau) \}}{\epsilon ^2}

como resultado de usar la inequidad de Markov

tiene maximo en τ=0

si se usa la inequidad de Cauchy-Schwarz:

E[XY]^2 \leq E[X^2]E[Y^2] R_X(\tau )^2 = E[X(t+ \tau) X(t)]^2 \leq E[X^2(t+ \tau)] E[X^2(t)] = R_X(0) |R_X(\tau)| \leq R_X(0)

Periódica en media cuadrática

si R_X(0) = R_X(d) , entonces R_X(\tau) es periódica con periodo d y X(t) es «periódica en media cuadrática».

R_X(\tau + d)| = R_X(\tau)

se aproxima a el cuadrado de la media cuando τ tiende a infinito

Proceso aleatorio estacionario

Referencia: León-García 9.6, 9.6.1 p518,521; Gubner 10.3, 10.2 p395, p392; Ross 10.7 p654, p656

Estacionario

Un proceso aleatorio discreto o contínuo en el tiempo 
X(t) es estacionario si la distribución conjunta de 
cualquier grupo de muestras No depende de la ubicación 
del tiempo de origen.
F_{X(t_1),..., X(t_k)} (x_1, ... x_k) = F_{X(t_1+\tau),..., X(t_k)} (x_1, ... x_k+\tau)

para cualquier τ, k, o muestras t1, … , tk.

Estacionario con cdf de 1er Orden

Un proceso aleatorio estacionario con cdf de primer orden debe ser independiente del tiempo

F_{X(t)} (x) = F_{X(t+\tau)} (x) = F_X(x)

para todo t, τ, dicho de otra forma, sus resultados son constantes

m_{X}(t) = E[X(t)]=m VAR[X(t)] = E[(X(t)-m)^2] = \sigma ^2

Estacionario con cdf de 2do Orden

Un proceso aleatorio estacionario con cdf de segundo orden debe depender solo de la diferencia de tiempo entre muestras y NO en particular del tiempo de las muestras

F_{X(t_1),X(t_2)} (x_1,x_2) = F_{X(0),X(t_2-t_1)} (x_1,x_2)

para todo t1, t2

R_{X}(t_1,t_2) = R_{X} (t_2-t_1) C_{X}(t_1,t_2) = C_{X} (t_2-t_1)

Estacionario en el Sentido Amplio (WSS o débil)

 Un proceso aleatorio discreto o contínuo en el tiempo 
X(t) es estacionario en el sentido amplio (WSS) 
si satiface que: 
m_X(t) = m C_{X}(t_1,t_2) = C_{X} (t_2-t_1)

De forma semejante, los procesos X(t) y Y(t) son estacionarios conjuntos si ambos son estacionarios en el sentido amplio y si covarianza cruzada depende solo de t1-t2.

R_{X}(t_1,t_2) = R_{X} (\tau) C_{X}(t_1,t_2) = C_{X} (\tau)

Autocovarianza PM

Referencia: León García Ejemplo 9.10 p495, Gubner Ejemplo 10.8 p389, Gubner Ejemplo 10.17 p396

Sea X(t) = cos(ω t + Φ), donde Φ es uniforme en el intervalo (-π,π) Encontrar la autocovarianza de X(t).

la variable aleatoria Φ tiene distribución uniforme en el intervalo, por lo que la función fΦ(φ) es constante = 1/[π – (-π)] = 1/2π.

Recordamos que:

E[g(x)] = \int_{-\infty}^{\infty} g(x) f(x) dx

Media (León-García 4.15 p158):

m_X(t) = E[cos(\omega t + \Phi)] = = \int_{-\pi}^{\pi} cos(\omega t + \Phi) \frac{1}{2\pi} d\Phi = \left. \frac{-1}{2\pi} (sin (\omega t + \Phi)) \right|_{-\pi}^{\pi} = \frac{-1}{2\pi} [sin (\omega t + (-\pi)) - sin (\omega t + \pi)] = 0

Autocovarianza

dado que el valor esperado es cero, la autocovarianza es igual a la autocorrelación

C_{X} (t_1,t_2) = R_X (t_1,t_2) - E[X(t_1,t_2)] = R_X (t_1,t_2) = E[cos(\omega t_1 + \Phi) cos(\omega t_2 +\Phi)]

Recordando que:

E[g(x)] = \int_{-\infty}^{\infty} g(x) f(x) dx cos(x) cos(y) = \frac{cos(x-y) + cos(x+y) }{2}

se tiene que:

= \int_{-\pi}^{\pi} [cos(\omega t_1 + \Phi) cos(\omega t_2 +\Phi)] \frac{1}{2\pi} d\Phi = \int_{-\pi}^{\pi} \frac{cos(\omega (t_1 - t_2))+cos(\omega (t_1 + t_2)+ 2\Phi)}{2} \frac{1}{2\pi} d\Phi = \int_{-\pi}^{\pi} \frac{cos(\omega (t_1 - t_2))}{2} \frac{1}{2\pi} d\Phi + \int_{-\pi}^{\pi} \frac{cos(\omega (t_1 + t_2 )+ 2\Phi)}{2} \frac{1}{2\pi} d\Phi

El primer integral, el coseno no depende de Φ, mientras que el segundo integral es semejante al intergral de la media y cuyo resultado es cero.

= \left. \frac{cos(\omega (t_1 - t_2))}{2} \frac{\Phi}{2\pi} \right|_{-\pi}^{\pi} + 0
C_{X} (t_1,t_2) = R_X (t_1,t_2) =
= \frac{1}{2} cos(\omega (t_1 - t_2))

Autocovarianza en AM

Referencia: León García Ejemplo 9.9 p495, Gubner Ej10.35 p496

Sea X(t) = A cos(2πt), donde A es una variable aleatoria, con un comportamiento semejante a la figura.

Encontrar el valor esperado , la autocorrelación y autocovarianza de de X(t).


El valor esperado se calcula a continuación, note que la media varia respecto a t y que el valor es cero para valores de t donde cos(2πt) =0.

E[X(t)] = E[A cos(2\pi t)] E[X(t)] = E[A] cos(2\pi t)

La autocorrelación es:

R_X(t_1,t_2) = E[A cos(2\pi t_1) A cos(2\pi t_2)] = E[A^{2} cos(2\pi t_1) cos(2\pi t_2)] = E[A^{2}] cos(2\pi t_1) cos(2\pi t_2)

usando:

2 cos(x)cos(y) = cos(x-y) + cos(x+y) cos(x)cos(y) = \frac{ cos(x-y) + cos(x + y)}{2}

se reemplaza:

= E[A^{2}] \frac{1}{2}[cos(2\pi t_1 - 2\pi t_2) + cos(2\pi t_1 + 2\pi t_2)] R_X(t_1,t_2) = E[A^{2}] \frac{[cos(2\pi (t_1 - t_2)) + cos(2\pi (t_1 + t_2))]}{2}

se observa que el valor de autocorrelación depende de las diferencias de tiempo t1 y t2.


La autocovarianza es:

Cov_X(t_1,t_2) = R_X(t_1,t_2) - E[X(t_1)]E[X(t_2)] = E[A^{2}] cos(2\pi t_1) cos(2\pi t_2) - E[A] cos(2\pi t_1)E[A] cos(2\pi t_2) = E[A^{2}] cos(2\pi t_1) cos(2\pi t_2) - E[A]^2 cos(2\pi t_1)cos(2\pi t_2) = (E[A^{2}] - E[A]^2) cos(2\pi t_1)cos(2\pi t_2) = Var[A] cos(2\pi t_1)cos(2\pi t_2)

con el mismo procedimiento de cos(x)cos(y):

Cov_X(t_1,t_2) = Var[A] \frac{[cos(2\pi (t_1 - t_2)) + cos(2\pi (t_1 + t_2))]}{2}

Autocorrelación, Autocovarianza con variable tiempo

Referencia: León-García 9.2.2 p493

Se puede usar los momentos de muestras en el tiempo para parcialmente especificar un proceso aleatorio al resumir la información contenida en las cdf conjuntas.


Procesos aleatorios Contínuos en tiempo

Media:

m_X(t) = E[X(t)] = \int_{-\infty}^{\infty} x f_{X(t)}(x) dx

Varianza:

VAR[X(t)] = \int_{-\infty}^{\infty} ( x - m_X(t))^2 f_{X(t)}(x) dx

donde f_{X(t)}(x) es la pdf de X(t). Note que ambas son funciones determinísticas de tiempo.

Autocorrelación

R_X(t_1,t_2) = E[X(t_1,t_2)] = = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} xy f_{X(t_1),X(t_2)}(x,y) dx dy

Autocovarianza:

C_X(t_1,t_2) = E[\{X(t_1) - m_X(t_1) \} \{ X(t_2) - m_X(t_2)\} C_X(t_1,t_2) = R_X(t_1,t_2) - m_X(t_1) m_X(t_2)

coeficiente de correlación

\rho (t_1, t_2) = \frac{C_x(t_1,t_2)}{\sqrt{C_x(t_1,t_1)}\sqrt{C_x(t_2,t_2)}}

El coeficiente de correlación es la medida en la cual una variable aleatoria puede predecirse como una función lineal de otra.


Procesos aleatorios Discretos en tiempo

Media:

m_X(n) = E[X(n)]

Varianza:

VAR[X(n)] = E[(X(n) - m_X(n))^2]

Autocorrelación

R_X(n_1,n_2) = E[X(n_1,n_2)]

Autocovarianza:

C_X(n_1,n_2) = E[\{X(n_1) - m_X(n_1) \} \{ X(n_2) - m_X(n_2)\} C_X(n_1,n_2) = R_X(n_1,n_2) - m_X(n_1) m_X(n_2)

 

Correlación

Referencia: León-García 5.6.2 p258, Gubner 2.4 p91

Correlación

La correlación entre dos variables aleatorias X y Y se define como E[XY].

La correlación determina cuando dos variables se encuentran linealmente relacionadas; es decir cuando una es función lineal de la otra.

R(X,Y) = E[XY]

Propiedades de la función de correlación

Simetría

R(X,Y) = R[Y,X] R(X,X) = E[X^2] \geq 0

Desigualdad de Cauchy-Schwarz

|R(X,Y)| = \sqrt{(E[X]^2 E[Y]^2)}

Covarianza

Retomando la función de covarianza de un proceso estocástico se muestra que:

Cov(X,Y) = E[XY] - E[X]E[Y] = = R(X,Y) - E[X]E[Y]

Coeficiente de correlación lineal

Al multiplicar una de las variables X o Y por un número se incrementa la covarianza, para una mejor medida se normaliza la covarianza y así tener los valores en una escala absoluta.

El coeficiente de correlación de X y Y se define por:

\rho_{X,Y} = \frac{Cov(X,Y)}{\sigma_X \sigma_Y} = \frac{E[XY]-E[X]E[Y]}{\sigma_X \sigma_Y} donde: \sigma_X=\sqrt{Var(X)} \sigma_Y=\sqrt{Var(Y)} -1 \leq \rho_{X,Y} \leq 1

El coeficiente de correlación es como máximo en magnitud 1.

Note que la correlación y el coeficiente de correlación no son lo mismo al resultar de la formula de covarianza.

Covarianza bivariada Ejemplo Ross 2_33

Referencia: Ejemplo Ross 2.33 p51

La función densidad conjunta de X,Y es

f(x,y) = \frac{1}{y} e^{-(y+x/y)} donde: 0<x , y<\infty

los intervalos también se expresan como: 0<x<∞  , 0<y<∞

a) Verifique que es una función de densidad conjunta

b) Encuentre la Cov(X,Y)


Solución propuesa

Para mostrar que f(x,y) es una función de densidad conjunta se debe mostar que es no negativa y que la integral doble es 1.

\int_{-\infty}^{\infty} \int_{-\infty}^{\infty}f(x,y) dy dx = = \int_{0}^{\infty} \int_{0}^{\infty} \frac{1}{y} e^{-(y+x/y)} dy dx = = \int_{0}^{\infty} e^{-y} \int_{0}^{\infty} \frac{1}{y} e^{-x/y} dx dy

Resolviendo parcialmente el integral

\int_{0}^{\infty} \frac{1}{y} e^{-x/y} dx = \frac{1}{y} \int_{0}^{\infty} e^{-x/y} dx = = \frac{1}{y} \left. \frac{e^{-x/y}}{-\frac{1}{y}}\right|_{0}^{\infty} = -[e^{-\infty /y} - e^{0} ] = -(0-1) = 1

con el resultado se continua con el integral anterior:

= \int_{0}^{\infty} e^{-y} dy = \left. \frac{e^{-y}}{-1} \right|_{0}^{\infty} = - (e^{-\infty} - e^{0}) = = -(0-1) = 1

con lo que se comprueba que es una función densidad de probabilidad conjunta.


Covarianza (X,Y)

Se obtendrán primero las funciones marginales para determinar el valor esperado de cada una de ellas:

f_Y(y) = \int_{-\infty}^{\infty} f(x,y) dx = \int_{0}^{\infty} \frac{1}{y} e^{-(y+x/y)} dx = = \int_{0}^{\infty} \frac{1}{y} e^{-y} e^{-x/y} dx = e^{-y} \int_{0}^{\infty} \frac{1}{y}e^{-x/y} dx = e^{-y}

usando parte del integral anterior, se conoce que es igual a 1

se requiere el valor esperado E[y]:

E[Y] = \int_{-\infty}^{\infty} y f(y) dy = \int_{0}^{\infty} y e^{-y} dy =

integración por partes dv=e^{-y}dy , v = - e^{-y} , u=y, du=dy

\int u dv = uv - \int v du E[Y] = \left. - ye^{-y} \right|_{0}^{\infty} - \int_{0}^{\infty} (-e^{-y}) dy = -(\infty e^{-\infty} - 0 e^{-0}) - \left. e^{-y} \right|_{0}^{\infty} = = 0 -[e^{-\infty}- e^{0}] = -[0 -1] E[Y] =1

El valor esperado E[x] se obtiene como:

E[x] = \int_{-\infty}^{\infty} \int_{-\infty}^{\infty} x f(x,y) dy dx = = \int_{0}^{\infty} \int_{0}^{\infty} x \frac{1}{y} e^{-(y+x/y)} dy dx = = \int_{0}^{\infty} e^{-y} \int_{0}^{\infty} \frac{x}{y} e^{-x/y} dx dy =

para el integral respecto a x:

\int_{0}^{\infty} \frac{x}{y} e^{-x/y} dx =

integración por partes dv=\frac{1}{y}e^{-x/y} dx , v = \frac{\frac{1}{y}}{-(1/y)} e^{-x/y}, u=x, du=dx

\int u dv = uv - \int v du = = \left. x(-e^{-x/y})\right|_{0}^{\infty} - \int_{0}^{\infty} - e^{-x/y} dx = = [- \infty ye^{-\infty/y} - (-0e^{0/y})] + \int_{0}^{\infty} e^{-x/y} dx = = [-0+0] + \left. \frac{e^{-x/y}}{-(1/y)}\right|_{0}^{\infty} = =-y[e^{-\infty/y} - e^{-0/y}] = y

reemplazando en E[x] y usando el resultado del integral de E[y]:

E[x] = \int_{0}^{\infty} y e^{-y} dy = 1

Para calcular el valor de E[XY]:

E[XY] = =\int_{-\infty}^{\infty} \int_{-\infty}^{\infty} xy f(x,y) dy dx = =\int_{0}^{\infty} \int_{0}^{\infty} xy \frac{1}{y} e^{-(y+x/y)} dy dx = = \int_{0}^{\infty} y e^{-y} \int_{0}^{\infty} \frac{x}{y} e^{-x/y} dx dy = = \int_{0}^{\infty} y e^{-y} (y) dy = \int_{0}^{\infty} y^2 e^{-y} dy =

nuevamente por partes: dv = e^{-y} dy, u=y^2 se obtiene

E[XY] = \int_{0}^{\infty} y^2e^{-y} dy = =\left. -y^2 e^{-y} \right|_{0}^{\infty} + \int_{0}^{\infty} 2ye^{-y}dy = 2E[Y] =2

En consecuencia:

Cov(X,Y) = E[XY] - E[X]E[Y] = 2-(1)(1) = 1

Covarianza bivariada

Referencia: Ross 2.5.3 p50

Covarianza

La covarianza y varianza de dos variables aleatorias X y Y se definen como:

Cov(X,Y) = E[(X-E[X])(Y- E[Y])] = E[XY] - E[X]E[Y]

Si X y Y son independientes, entonces Cov(X,Y) = 0

Propiedades

  1. Cov(X,X) = Var(X)
  2. Cov(X,Y) = Cov(Y,X)
  3. Cov(cX,Y) = c Cov(X,Y)
  4. Cov(X,Y+Z) = Cov(X,Y) + Cov(X,Z)