In [1]:
voir le code
Les séries temporelles ont souvent une composante périodique, comme par exemple la variation journalière de température.
Read data files¶
In [4]:
# Read
data = np.load('./timeseriesdata.npy').item()
print data.keys()
In [9]:
#key = './data/Tfrigo.csv'
key = './data/ICNSA.csv'
X, Y = data[key]
plt.figure(figsize=(14, 4) )
plt.plot( Y, 'k-', alpha=0.5 );
plt.xlabel( 'X' ); plt.title(key); plt.xlim([1, len(Y)]);
Autocorrelation¶
In [10]:
def autocorrelation(x):
"""
Compute the autocorrelation of the signal, based on the properties of the
power spectral density of the signal.
http://stackoverflow.com/a/40154897
"""
xp = x - np.mean(x)
f = np.fft.fft(xp)
p = np.array([np.real(v)**2+np.imag(v)**2 for v in f])
p = np.absolute( f )**2
pi = np.fft.ifft(p)
return np.real(pi)[:x.size/2]/np.sum(xp**2)
On cherche les changements de signe (de positif vers negatif) de la dérivée pour identitfier les maximums :
In [11]:
def findargmax( dY ):
signe = np.sign( dY[0] )
argmax = []
for i, dy in enumerate( dY ):
new_sign = np.sign( dy )
if signe == +1 and new_sign==-1:
argmax.append( i )
signe = new_sign
return argmax
Test de repliement de periode¶
In [12]:
autocor = autocorrelation( Y )
xMax = 250
# find maximums
argmax = findargmax( np.diff( autocor )[:xMax] )
maxs = [ autocor[i] for i in argmax ]
print argmax
plt.figure(figsize=(14, 7) )
plt.subplot( 2, 1, 1 )
plt.plot( np.diff( autocor[:xMax] ) , '-r' );
plt.axhline(y=0., color='k', linestyle='-');
plt.subplot( 2, 1, 2 );
plt.plot( autocor[:xMax] );
plt.plot( argmax, maxs, 'or' );
La periode principal est $52$. Elle correspond bien au 52 semaines d'une année. La seconde seconde période étant une demi année.
In [8]:
period = 52
offset = 0
n = int( float(len(Y)-offset)/float(period) )
smoothX, smoothY = [], []
smooth = np.zeros( np.size(Y))
plt.figure(figsize=(14, 4) )
for i in range(n):
Yi = Y[ offset+i*period : offset+(i+1)*period ]
Yi_zero = Yi - np.mean( Yi )
plt.plot( Yi_zero, 'k-', alpha=0.3 );
smoothX.append( (offset+i*period +offset+(i+1)*period)/2 )
smoothY.append( np.median( Yi ) )
smooth[ offset+i*period : offset+(i+1)*period ] = Y[ offset+i*period : offset+(i+1)*period ] - np.median( Yi )
plt.xlim([0, period-1]);
plt.figure(figsize=(14, 4) )
plt.plot( smoothX, smoothY, 'r' )
plt.plot( Y, 'k', alpha=0.3 );
plt.xlim([0, len(Y)-1]);
#plt.figure(figsize=(14, 4) )
#plt.plot( smooth, 'b' )
Out[8]:
