%% Proceso estocástico: media, varianza, correlación temporal
% © 2019 Antonio Sala Piqueras, Universitat Politècnica de València. Todos los derechos reservados.
% Presentación en vídeo en http://personales.upv.es/asala/YT/V/estoctm.html
% Este código ejecutó sin errores en Matlab R2018b (Linux)
% Tabla de Contenidos
% Generación de datos
% Análisis de los datos
%     Media, varianza y desviación típica
%     Covarianza entre distintos instantes
%     Correlación entre distintos instantes
% Funciones auxiliares
%
%% Generación de datos

MuestrasExperimento=400;
datos=formulasecreta(MuestrasExperimento);
plot(datos), grid on, title('Una realización'), xlabel('tiempo'), ylabel('Volt')

NExperimentos=2000;
MuchosDatos=zeros(NExperimentos,MuestrasExperimento);
for i=1:NExperimentos
    MuchosDatos(i,:)=formulasecreta(MuestrasExperimento);
end
plot(MuchosDatos'), grid on, title('2000 realizaciones'), xlabel('tiempo'), ylabel('Volt')

%% Análisis de los datos
%% Media, varianza y desviación típica

Media=mean(MuchosDatos);
plot(Media), grid on, title('MEDIA'), xlabel('tiempo'), ylabel('Volt')
Varianza=var(MuchosDatos);
DsvT=sqrt(Varianza);
plot(DsvT), grid on, title('Desv. Típica'), xlabel('tiempo'), ylabel('Volt')

%% Covarianza entre distintos instantes

MatrizVCCompleta=cov(MuchosDatos);
plot(MatrizVCCompleta(1,:)), grid on, xlabel('tiempo'), title('Covarianza entre x(1) y resto.'), ylabel('Volt^2')
plot(MatrizVCCompleta(20,:)), grid on, xlabel('tiempo'), title('Covarianza entre x(20) y resto.'), ylabel('Volt^2')
plot(MatrizVCCompleta(50,:)), grid on, xlabel('tiempo'), title('Covarianza entre x(50) y resto.'), ylabel('Volt^2')
plot(MatrizVCCompleta(250,:)), grid on, xlabel('tiempo'), title('Covarianza entre x(250) y resto.'), ylabel('Volt^2')
plot(MatrizVCCompleta(350,:)), grid on, xlabel('tiempo'), title('Covarianza entre x(350) y resto.'), ylabel('Volt^2')
TT=[5 10 30 50 100 150 200 300 400];
plot(MatrizVCCompleta(TT,:)'), grid on, xlabel('tiempo'), ylabel('Volt^2')

%% Correlación entre distintos instantes

MatrCORR=inv(diag(DsvT))*MatrizVCCompleta*inv(diag(DsvT));
plot(MatrCORR(TT,:)'), grid on, xlabel('tiempo'), title('Correlación entre algunos instantes y resto.'), ylabel('\rho (adimensional)')







%
%% Funciones auxiliares

function D=formulasecreta(MuestrasExp)
    x=10+0.1*randn();
    D=zeros(1,MuestrasExp);
    for i=1:MuestrasExp
        D(i)=x;
        x=0.978*x+0.022*randn()*11;
    end
   
end