
% PROGRAMMA ScatterPlot
% Autore: Claudio Prati
% Politecnico di Milano, 28 Aprile 2025

% Data una serie temporale generata dal modello x(n+1)=rho*x(n)+kk*w(n)
% dove rho è il coefficiente di correlazione a passo 1 del processo casuale, kk è un 
% coefficiente che garantisce la stazionarietà del processo casuale e w(n) è un processo casuale
% indipendente da x(n) e con ddp gaussiana, si visualizza lo scatterplot
% x(n) X(n+1).

% Viene inoltre visualizzata tramite istogramma 3D la densità di
% probabilità congiunta tra x(n+1) e x(n).
% Da questa è facile realizzare l'andamento della ddp di x(n+1)
% condizionata a x(n) e il suo valore medio che corrisponde alla predizione
% che minimizza l'errore quadratico medio tra dato previsto e dato vero.


clear all
close all

ro=input('coefficiente di correlazione? ');
N=2e7;
M=500;
% Calcolo di kk per garantire la stazionarietà del processo
kk=sqrt(1-ro^2);
x=zeros(1,N+1);
w=randn(1,N);
x(1)=randn(1,1);

% Generazione del processo casuale
for n=1:N;
    x(n+1)=ro*x(n)+kk*w(n);
end


% Visualizzazione dello scatterplot di M campioni 

figure(1)
for n=1:M
plot(x(n),x(n+1),'.')
axis('image')
pause(.001)
hold on
end

hold off
pause




% visualizzazione della densità di propabilità del processo x(n)
figure(2)
da=0.05;
cnt=hist(x(1:N),[-5:da:5]);
ddpx=cnt/N/da;
C=(-5:da:5);
plot(C,ddpx)
title('Ddp del processo casuale x(n)')

pause
% visualizzazione della densità di probabilità congiunta tra x(n) e x(n+1)
figure(3)
yx(1:N,2)=x(2:N+1);
yx(1:N,1)=x(1:N);
da=.05;
db=.05;
cnt=hist3(yx,{-5:da:5 -5:db:5});
ddpyx=cnt/N/(da*db);
imagesc((-5:da:5),(-5:db:5),ddpyx)
axis('image')
colorbar

title('Ddp congiunta di x(n+1) e x(n)')


xn=0;

% visualizzazione della densità di propabilità di x(n+1) condizionata a x(n)
for k=1:5;
xn=input(' Inserire un valore di x(n) compreso tra -3 e +3 ... >>  ');
nn=floor((xn+5)/da);
figure(4)
plot((-5:db:5),ddpyx(nn,:)/(ddpx(nn)))
grid on
pause
end
