DHOB (IU5SGN): gennaio 2024

venerdì 12 gennaio 2024

Rimozione trend da serie tempo

Una prova per rimuovere la componente stagionale da una serie tempo usando la libreria StatsModels di Python

In alto dati originali. al centro trend, seguend dati periodici stimati ed in basso calcolo dei residui

#!/usr/bin/env python
# coding: utf-8

# In[1]:

from IPython.display import display

import numpy as np
import pandas as pd
pd.set_option('display.max_rows', 15)
pd.set_option('display.max_columns', 500)
pd.set_option('display.width', 1000)

import matplotlib.pyplot as plt
from datetime import datetime
from datetime import timedelta
from pandas.plotting import register_matplotlib_converters
from mpl_toolkits.mplot3d import Axes3D

from statsmodels.tsa.stattools import acf, pacf
from statsmodels.tsa.statespace.sarimax import SARIMAX
import statsmodels.api as sm
import statsmodels.formula.api as smf

register_matplotlib_converters()
from time import time
import seaborn as sns
sns.set(style="whitegrid")

from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.cluster import KMeans
from sklearn.covariance import EllipticEnvelope

import warnings
warnings.filterwarnings('ignore')

RANDOM_SEED = np.random.seed(0)

# # Carica i dati

# In[2]:

def parser(s):
    return datetime.strptime(s, '%Y-%m-%d %H:%M:%S')

# In[3]:

#read data
est = pd.read_csv('/home/luca/luca/liscione.csv',sep=";",parse_dates=[0], index_col=0, date_parser=parser,header=0)
print(est)

# In[4]:

est = est.asfreq(freq='1H', method='ffill')

# ### Introduce an Anomaly

# In[5]:

idx = pd.IndexSlice

# In[6]:

start_date = datetime(2023,10,5,0,0,0)
end_date = datetime(2023,12,31,23,59,0)
lim_est = est[start_date:end_date]

# In[7]:

lim_est

# In[8]:

del lim_est['Temp']

# In[9]:

plt.figure(figsize=(10,4))
plt.plot(lim_est)
plt.title('Est', fontsize=20)
plt.ylabel('Est(mm)', fontsize=16)
for year in range(start_date.year,end_date.year):
    plt.axvline(pd.to_datetime(str(year)+'-01-01'), color='k', linestyle='--', alpha=0.2)

# ## Seasonal Decompose

# In[10]:

from statsmodels.tsa.seasonal import seasonal_decompose
import matplotlib.dates as mdates

# In[18]:

plt.rc('figure',figsize=(12,8))
plt.rc('font',size=15)

result = seasonal_decompose(lim_est,model='additive')
print(result.trend)
plt.savefig('trend')

fig = result.plot()
fig.savefig('trend.png')   # save the figure to file
result.trend.to_csv('/home/luca/tempo.csv', date_format='%Y-%m-%d_%H',sep=";")

# In[12]:

plt.rc('figure',figsize=(12,6))
plt.rc('font',size=15)

fig, ax = plt.subplots()
x = result.resid.index
y = result.resid.values
print(x.size)
ax.plot_date(x, y, color='black',linestyle='--')
ax.set_ylim([-0.7, 0.7])

fig.autofmt_xdate()
plt.savefig('residual')

plt.show()

# In[13]:

temp = est[start_date:end_date]

del temp['Est']
print(temp.index)
temp.index.to_csv('tempo2.csv', date_format='%Y-%m-%d_%H',delimiter=";")

#np.savetxt('tempo.csv',temp.index,delimiter=";")

# In[ ]:

# In[ ]:

fig, ax1 = plt.subplots()

color = 'tab:red'
ax1.set_xlabel('time')
ax1.set_ylabel('Temp', color=color)
ax1.plot(temp, color=color)
ax1.tick_params(axis='y', labelcolor=color)

ax2 = ax1.twinx()  # instantiate a second axes that shares the same x-axis

color = 'tab:blue'
ax2.set_ylabel('Est', color=color)  # we already handled the x-label with ax1
ax2.plot(result.trend, color=color)
ax2.tick_params(axis='y', labelcolor=color)

plt.savefig('sovrapposto')

plt.show()

# In[ ]:

np.savetxt('trend.csv',result.trend.values,delimiter=";")
np.savetxt('temp.csv',temp.values,delimiter=";")

# In[ ]:

import matplotlib.animation as animation
fig,ax = plt.subplots()

scat = ax.scatter(result.trend.values, temp.values)
plt.gca().invert_xaxis()
ax.set(xlabel='Estensimetro (mm)')
ax.set(ylabel='Temperatura (C)')

def update(i):
    scat.set_offsets((result.trend.values[i],temp.values[i]))
    return scat,

ani = animation.FuncAnimation(fig,update, repeat=True, frames=len(temp)-1, interval=10)
#writer = animation.PillowWriter(fps=15)
#ani.save('animation.gif',writer=writer)
#plt.savefig('scatterplot')

plt.show()

Animazione di grafico con Matplotlib

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
import pandas as pd
df = pd.read_csv("tempo.csv", dtype=str,sep=";")
print(df)
my_data = np.genfromtxt('./scatterplot.csv', delimiter=';')
x = my_data[:,1]
y = my_data[:,0]
fig = plt.figure()
ax = fig.add_subplot()
plt.ylim(-5, 22)
plt.ylabel ('Temperatura (C)')
plt.xlim(-59, -53)
plt.xlabel ('Estensimetro (mm)')
tempo = ax.text(0.8, 0.9, '', horizontalalignment='center',verticalalignment='center', transform=ax.transAxes)
graph, = plt.plot([], [], 'o', markersize=1)
def animate(i):
graph.set_data(x[:i+1], y[:i+1])
tempo.set_text(df.Data[i])
return graph,
ani = FuncAnimation(fig, animate, frames=2000, interval=50)
#FFwriter = animation.FFMpegWriter(fps=10)
#ani.save('animation.mp4')
plt.show()