[PYTHON] Test de Dicky Fuller par des modèles statistiques
Cet article a été retiré de ce que j'avais publié il y a quelques jours et a été amélioré à nouveau.
Il est important de savoir si les données de séries chronologiques à analyser suivent une marche aléatoire à partir des deux points suivants.
- La tendance qui s'est produite est-elle une fausse tendance ou la tendance se produit-elle?
- Lorsque les deux séries chronologiques montrent une corrélation, s'agit-il d'une fausse corrélation ou se produisent-elles? C'est le but.
De nombreuses séries chronologiques économiques et niveau de stock (série originale) suivent une marche aléatoire. Une série chronologique qui suit une marche aléatoire est une série chronologique non stationnaire dans laquelle la moyenne et la variance diffèrent des précédentes au fil du temps. Dans une marche aléatoire, la moyenne change de manière probabiliste avec le temps et la variance augmente proportionnellement à la racine carrée du passage du temps. Il existe également diverses séries chronologiques non stationnaires. L'un d'eux est la stabilité de la tendance dans laquelle la valeur moyenne augmente et diminue régulièrement avec le temps. Pour déterminer s'il s'agit d'une marche aléatoire, déterminez si la série chronologique a une racine unitaire. Il existe un test Dicky Fuller comme méthode. Cependant, il est souvent souligné que la puissance de détection est faible. Je vais donc l'examiner un peu.
Avant de passer au test DF
Tout d'abord, un modèle d'auto-retour du premier ordre (AR (1))
Puisque la série temporelle diverge dans $ \ rho> 1 $, définissez $ \ rho <1 $. Ensuite, le modèle AR (1) fera de même.
Vérifions jusqu'à ce point en utilisant python.
import numpy as np
import matplotlib.pyplot as plt
def generate(y0,rho,sigma,nsample):
e = np.random.normal(size=nsample)
y=[]
y.append(rho*y0+e[0])
for i,u in enumerate(e[1:]):
y.append(rho*y[i]+u)
plt.plot(y)
generate(1,1,1,100)
Ensuite, définissons $ rho = 0,9 $.
generate(1,0.9,1,100)
Vous pouvez voir qu'il converge vers une certaine valeur. Ensuite, définissons sigma = 0 pour les deux.
generate(1,1,0,100)
generate(1,0.9,0,100)
Il reste 1 sur la marche aléatoire et converge vers zéro sur AR (1). Ceci n'est pas pratique, alors réécrivez le modèle AR (1) comme suit.
En supposant que la valeur moyenne de $ y_t $ est $ \ mu $, si n est assez grand, alors $ \ sum_ {t = 2} ^ nu_t = 0 $
Cela devient $ \ mu = c + \ rho \ mu $. Donc
def generate(y0,rho,c,sigma,nsample):
e = np.random.normal(size=nsample)
y=[]
y.append(rho*y0+e[0]*sigma+c)
for i,u in enumerate(e[1:]):
y.append(rho*y[i]+u*sigma+c)
plt.plot(y)
generate(1,0.9,0.1,0,100)
est devenu. La valeur moyenne de AR (1) est désormais constante dans le temps. Rendons sigma plus petit et vérifions-le.
generate(1,0.9,0.1,0.01,100)
D'une manière ou d'une autre, vous pouvez voir qu'il revient au centre. Le modèle AR (1) devient stable (faible stable) si certaines conditions sont remplies. Ce c est parfois appelé le taux de dérive. La vitesse de dérive est la vitesse à laquelle la moyenne change. Mais ça ne change pas ici. Cependant, c est Si c'est $ c> \ mu (1- \ rho) $ ou $ c <\ mu (1- \ rho) $, le taux changera.
Dans la théorie stochastique, la dérive stochastique est le changement de la valeur moyenne du processus stochastique. Il existe un taux de dérive avec un concept similaire, mais cela ne change pas de manière probabiliste. C'est un paramètre non probabiliste.
Dans les études à long terme sur les événements à long terme, les caractéristiques des séries chronologiques sont souvent conceptualisées en les décomposant en composantes tendance, périodiques et probabilistes (dérives probabilistes). Les composantes de la dérive périodique et probabiliste sont identifiées par l'analyse d'autocorrélation et la différence des tendances. L'analyse d'autocorrélation permet d'identifier la phase correcte du modèle ajusté et la répétition de la différence transforme la composante de dérive stochastique en bruit blanc.
Tendance stable
Le processus de tendance stable $ y_t $ suit $ y_t = f (t) + u_t $. Où t est le temps, f est une fonction déterministe et $ u_t $ est une variable de probabilité stationnaire avec une moyenne à long terme de zéro. Dans ce cas, le terme stochastique est stationnaire et il n'y a donc pas de dérive stochastique. La composante déterministe f (t) n'a pas de valeur moyenne qui reste constante sur le long terme, mais la série temporelle elle-même peut dériver. Cette dérive déterministe peut être calculée en renvoyant $ y_t $ à t et en la supprimant des données une fois qu'un résidu stable est obtenu. Le processus de stabilisation de tendance le plus simple est $ y_t = d \ cdot t + u_t $. $ d $ est une constante.
Échelonné
Le processus de racine unitaire (stationnaire) suit $ y_t = y_ {t-1} + c + u_t $. $ u_t $ est une variable de probabilité stationnaire avec une valeur moyenne qui sera nulle à long terme. Où c est un paramètre de dérive non probabiliste. Même sans le choc aléatoire $ u_t $, la moyenne de y change de c dans le temps. Cette non-stationnarité peut être supprimée des données en prenant une différence de premier ordre. La variable graduée $ z_t = y_t-y_ {t-1} = c + u_t $ a une moyenne constante à long terme de c, il n'y a donc pas de dérive dans cette moyenne. Ensuite, considérons c = 0. Ce processus racine unitaire est $ y_t = y_ {t-1} + c + u_t $. Bien que $ u_t $ soit un processus stationnaire, la présence de son choc probabiliste $ u_t $ provoque une dérive, notamment probabiliste, dans $ y_t $. Une fois qu'une valeur u non nulle se produit, elle est incluse dans y pour la même période. Ce sera une période derrière y après une période et affectera la valeur y pour la période suivante. y affecte les nouvelles valeurs y les unes après les autres. Par conséquent, après que le premier choc affecte y, sa valeur est incorporée en permanence dans la moyenne de y, résultant en une dérive probabiliste. Donc, cette dérive peut également être éliminée en différenciant d'abord y pour obtenir z qui ne dérive pas. Cette dérive stochastique se produit également dans $ c \ ne0 $. Cela arrive donc aussi avec $ y_t = y_ {t-1} + c + u_t $. J'ai expliqué que c est un paramètre de dérive non probabiliste. Pourquoi n'est-ce pas une dérive définitive? En fait, la dérive qui résulte de $ c + u_t $ fait partie du processus stochastique.
numpy.random.normal De la référence numpy
numpy.random.normal(loc=0.0, scale=1.0, size=None)
loc float ou array_like of floats: valeur moyenne de la distribution
scalefloat ou array_like of floats: écart type de la distribution
est. Considérons maintenant la signification de loc. La valeur moyenne des échantillons générés à partir de cette distribution suit de près $ N (μ, σ ^ 2 / n) $. Voyons voir comment ça fonctionne. Osez définir loc = 1, échelle = 1.0. J'ai essayé de trouver la moyenne et l'écart type de l'échantillon avec n = 2. Cependant, il existe deux méthodes de calcul. L'un est Trouvez la moyenne et l'écart type de $ z_t $. L'autre est Soustrayez 1 de chaque observation $ z_t $ pour trouver la moyenne et l'écart type de l'échantillon. Ajoutez ensuite 1 à la moyenne de l'échantillon.
import numpy as np
import matplotlib.pyplot as plt
s = np.random.normal(1, 1, 2)
print(np.mean(s),np.std(s))
o=np.ones(len(s))
Print(np.mean(s-o)+1,np.std(s-o))
# 0.9255104221128653 0.5256849930003175
# 0.9255104221128653 0.5256849930003175 #1 a été soustrait de l'échantillon et 1 a été ajouté à la valeur moyenne.
Le résultat est le même. On peut voir que np.random.normal (1,1, n) génère des nombres aléatoires avec 1 / n + np.random.normal (0,1, n). Mais c'est une histoire dans l'ordinateur. En réalité, il devrait être impossible de distinguer celui qui crée la dérive stochastique. Générons donc des données 1000000 fois pour n = 2 et n = 100 et voyons à quoi ressemble le résultat.
Distribution des moyennes arithmétiques d'échantillons de processus stables échelonnés
La variable graduée $ z_t = y_t-y_ {t-1} = c + u_t $ a une moyenne constante à long terme de c. La moyenne de l'échantillon suit de près $ N (μ, σ ^ 2 / n) $ si n est grand. Par conséquent, la valeur moyenne n'est pas constante. Cependant, dans ce cas, il n'y a pas de dérive de la valeur moyenne. Il a une valeur moyenne constante sur le long terme.
m1=[]
trial=1000000
for i in range(trial):
s = np.random.normal(1, 1, 2)
m1.append(np.mean(s))
m2=[]
for i in range(trial):
s = np.random.normal(1, 1, 100)
m2.append(np.mean(s))
bins=np.linspace(np.min(s),np.max(s),1000)
plt.hist(m1,bins=bins,label='n=2')
plt.hist(m2,bins=bins,label='n=100',alpha=0.5)
plt.legend()
plt.show()
La moyenne de l'échantillon et l'écart type forment une distribution et ne sont pas déterminés de manière unique. N doit être beaucoup plus grand pour être unique. Par conséquent, les caractéristiques de la série chronologique sont décomposées en composantes de tendance, composantes périodiques et composantes probabilistes (dérive probabiliste) pour les conceptualiser, mais il faut faire attention lors de la gestion de la dérive probabiliste et du taux de dérive.
Une dérive probabiliste se produit dans le processus de racine unitaire. Même s'il existe un taux de dérive défini, il est impossible de le distinguer.
Distribution de y avec racines unitaires le dernier jour
Dans le processus de racine unitaire, l'écart type augmente à mesure que le nombre d'échantillons augmente. La valeur moyenne de y avec une taille d'échantillon de n est zéro à c = 0, mais la distribution augmente à mesure que l'écart-type $ \ sqrt {n \ cdot \ sigma ^ 2} $ augmente.
m1=[]
trial=1000000
for i in range(trial):
s = np.cumsum(np.random.normal(0, 1, 2))
m1.append(s[-1])
m2=[]
for i in range(trial):
s = np.cumsum(np.random.normal(0, 1, 100))
m2.append(s[-1])
bins=np.linspace(np.min(m2),np.max(m2),1000)
plt.hist(m1,bins=bins,label='n=2')
plt.hist(m2,bins=bins,label='n=100',alpha=0.5)
plt.legend()
plt.show()
La distribution est différente de la précédente. Plus la taille de l'échantillon est grande, plus la distribution est large. Essayons de tracer une série.
plt.plot(s)
Dans ce cas, le niveau de chaque échantillon change par étapes. Ce changement est dû à un choc stochastique $ u_t $. Ce à quoi il a été ajouté est tracé. Une série de dérives stochastiques est tracée. Cela peut former une grande tendance. Il correspond à la zone aux deux extrémités de la carte de répartition de n = 100. Quelle plage dépend de la définition de la tendance. Cette tendance est une tendance probabiliste. Cela se distingue clairement des tendances déterministes créées par la stabilité des tendances.
Processus de retour automatique de premier ordre
Le processus de retour automatique de premier ordre peut être réécrit comme suit.
Lorsque $ \ rho = 1 $, ce sera un processus de marche aléatoire.
#Initialisation
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
from scipy.stats import norm
from statsmodels.tsa.stattools import adfuller
#np.random.seed(0)
# AR(1)Génération de processus
#x0 est le prix initial
#nsample est le nombre de périodes dans un essai
# rho=1: random walk
## c=c1-rho est la dérive, c1=Lorsqu'il est 1, il n'y a pas de dérive Marche aléatoire d est une tendance définie
# rho<Modèle à retour automatique 1 au 1er ordre
##Les effets de x0, c et d disparaissent avec le temps
def ar1(rho,x0,c1,d,sigma,nsample):
e = np.random.normal(size=nsample)#random variable
rhoo=np.array([rho ** num for num in range(nsample)]) #effects of rho
rhoo_0=np.array([rho**(num+1) for num in range(nsample)]) #effects of rho
t=np.arange(nsample) #deterministic trend
one=np.ones(nsample) #effect of x0
x00=x0*one*rhoo_0 #effect of rho to x0
c=c1-rho # constant of autoregressive mdoel
c0=np.cumsum(c*one*rhoo) #effect of rho to c
d0=(d*t*rhoo) # effect of rho to determinist trend
y=np.cumsum(rhoo*e*sigma)+x00+c0+d0
return y
Distribution de y le dernier jour du processus AR (1)
m1=[]
trial=1000000
nsample=250
x0=1
for i in range(trial):
s = np.cumsum(np.random.normal(0, 1, nsample))
m1.append(s[-1]+x0)
def generate(rho,x0,c1,sigma,nsample,trial):
m=[]
for i in range(trial):
e = np.random.normal(0, 1, nsample)
rhoo=np.array([rho ** num for num in range(nsample)]) #effects of rho
rhoo_0=np.array([rho**(num+1) for num in range(nsample)]) #effects of rho
one=np.ones(nsample) #effect of x0
x00=x0*one*rhoo_0 #effect of rho to x0
c=c1-rho # constant of autoregressive mdoel
c0=np.cumsum(c*one*rhoo) #effect of rho to c
s=np.cumsum(rhoo*e*sigma)+x00+c0
m.append(s[-1])
return m
c1=1
x0=1
rho=0.99
sigma=1
m2=generate(rho,x0,c1,sigma,nsample,trial)
rho=0.9
m3=generate(rho,x0,c1,sigma,nsample,trial)
rho=0.5
m4=generate(rho,x0,c1,sigma,nsample,trial)
bins=np.linspace(np.min(m1),np.max(m1),1000)
plt.hist(m1,bins=bins,label='rw')
plt.hist(m2,bins=bins,label='ar1 rho=0.99',alpha=0.5)
plt.hist(m3,bins=bins,label='ar1 rho=0.9',alpha=0.5)
plt.hist(m4,bins=bins,label='ar1 rho=0.5',alpha=0.5)
plt.legend()
plt.show()
Dans le modèle de retour automatique du premier ordre, vous pouvez voir que lorsque $ \ rho $ devient inférieur à 1, il erre autour de la valeur initiale. A mesure que $ \ rho $ s'approche de 1, il devient difficile de le distinguer d'une marche aléatoire. La figure suivante est pour n = 7. Dans environ 7 jours, la propagation de tout $ \ rho $ sera similaire.
Détection du pouvoir (Introduction aux tests d'hypothèses statistiques par des modèles de statistiques)
Dans le test d'hypothèse, l'hypothèse nulle et l'hypothèse alternative sont faites. L'hypothèse nulle consiste en une population. Comme la population est généralement inconnue, une estimation peut être utilisée. Les spécimens sont traités comme extraits de cette population. Dans la figure suivante, la répartition de la population est indiquée par la ligne pointillée bleue. L'échantillon donne la distribution d'échantillon de cette statistique. Il est indiqué par une ligne orange. Comparez les deux. Si les deux distributions sont proches, il est difficile de juger que les deux distributions sont différentes, et si elles sont éloignées, on peut juger qu'elles sont différentes.
Ce jugement peut être fait en utilisant le niveau de signification. Le niveau de signification $ \ alpha $ est la zone la plus à droite de la population. 10%, 5%, etc. sont utilisés. La zone à droite de la ligne verticale bleue dans la distribution bleue. La ligne verticale bleue montre la limite où l'hypothèse nulle est rejetée. La zone à l'intérieur de cette ligne est appelée zone d'acceptation et la zone à l'extérieur de cette ligne est appelée zone de rejet.
Une fois l'échantillon obtenu, des statistiques telles que la moyenne et la variance peuvent être calculées. La valeur p est calculée à partir de la distribution de la population, qui est la probabilité d'être rejeté lorsque l'hypothèse nulle est rejetée avec la statistique comme limite inférieure. Ensuite, sur la figure, il s'agit de la zone à droite de la ligne verticale rouge. Si cette valeur est inférieure au critère de signification, on estime que l'échantillon (statistique) a été obtenu à partir d'une distribution différente de la population (population). Sur la figure, la ligne verticale rouge se trouve à gauche de la ligne verticale bleue, mais si la valeur moyenne de la distribution de l'échantillon orange est à droite de la ligne verticale bleue, il est clair que sa valeur p est inférieure au niveau de signification.
La capacité de rejeter correctement une fausse hypothèse nulle est appelée détection. La puissance de détection $ \ beta $ correspond à la zone à droite de la distribution orange des lignes verticales bleues. Si la partie $ \ beta $ est grande, la probabilité que l'hypothèse alternative soit correcte augmentera lorsque l'hypothèse nulle est rejetée.
Afin d'augmenter ce $ \ beta $, il est possible que les deux distributions soient éloignées, ou que les positions des centres des distributions soient significativement différentes. L'autre est la grande taille du spécimen. Ensuite, la forme de la distribution sera rétrécie, et si la taille de l'échantillon est suffisamment grande pour considérer la distribution comme une ligne, alors les deux lignes droites seront comparées.
Test de Dicky-Fuller
Un processus AR (1) simple
$ y_t= \rho y_{t-1}+u_t$
Et soustrayez $ y_ {t-1} $ des deux côtés
Il existe trois types.
- Marche aléatoire sans dérive
\Delta y_t=\delta y_{t-1}+u_t - Marche aléatoire avec dérive
\Delta y_t=c+\delta y_{t-1}+u_t - Marche aléatoire avec dérive + tendance temporelle
\Delta y_t=c+c\cdot t+\delta y_{t-1}+u_t Ce sont des différences de premier ordre sur le côté gauche, ce qui signifie que les différences ont des dérives ou des tendances temporelles.
# Dicky-Test plus complet
def adfcheck(a0,mu,sigma,nsample,trial):
prw=0
prwc=0
prwct=0
for i in range(trial):
y=ar1(a0,1,mu,sigma,nsample)
rw=sm.tsa.adfuller(y,maxlag=0,regression='nc')[1]
rwc=sm.tsa.adfuller(y,maxlag=0,regression='c')[1]
rwct=sm.tsa.adfuller(y,maxlag=0,regression='ct')[1]
if rw>0.05:
prw+=1
if rwc>0.05:
prwc+=1
if rwct>0.05:
prwct+=1
print('nsample',nsample,prw/trial,prwc/trial,prwct/trial)
Test de marche de Landau
Marche aléatoire sans dérive
rho=1
c1=1
c1=(c1-1)/250+1
d=0
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.9448 0.8878 0.8528
# nsample 20 0.9635 0.9309 0.921
# nsample 250 0.9562 0.9509 0.9485
Marche aléatoire avec dérive
rho=1
c1=1.1
c1=(c1-1)/250+1
d=0
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.9489 0.8889 0.8524
# nsample 20 0.9649 0.9312 0.921
# nsample 250 0.9727 0.9447 0.9454
Marche aléatoire + dérive + tendance déterministe
rho=1
c1=1.1
c1=(c1-1)/250+1
d=0.9
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.9971 0.948 0.8581
# nsample 20 1.0 1.0 0.983
# nsample 250 1.0 1.0 0.9998
Si la marche aléatoire doit être rejetée
Modèle AR (1) sans dérive
rho=0.99
c1=rho
#c1=(c1-1)/250+1
d=0
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.7797 0.9012 0.8573
# nsample 20 0.6188 0.9473 0.9225
# nsample 250 0.0144 0.7876 0.944
Modèle AR (1) avec dérive
rho=0.90
c1=1.1
c1=(c1-1)/250+1
d=0
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.9715 0.8868 0.8539
# nsample 20 0.9989 0.9123 0.912
# nsample 250 1.0 0.0306 0.1622
Modèle AR (1) + dérive + tendance déterministe
rho=0.90
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.999 0.9404 0.8576
# nsample 20 1.0 1.0 0.918
# nsample 250 1.0 1.0 0.0044
Nous pouvons voir que nous devons utiliser un modèle avec les spécifications correctes, que la taille de l'échantillon doit être grande et que ρ doit être éloigné de 1 lors du rejet de l'hypothèse nulle.
À propos de l'influence du bruit
Utilisons la distribution de Cauchy pour générer un gros bruit de queue et utilisons-la pour voir comment elle affecte la série chronologique. Regardons d'abord la distribution de la distribution de Cauchy.
from scipy.stats import cauchy
trial=1000000
m1=[]
for i in range(trial):
s = np.random.normal(0, 1, 250)
m1.append(np.mean(s))
m2=[]
for i in range(trial):
s=[]
i=1
while i<=250:
s0 = (cauchy.rvs(0,1,1))
if abs(s0)<10:
s.append(s0)
i+=1
m2.append(np.mean(s))
bins=np.linspace(np.min(m2),np.max(m2),1000)
plt.hist(m1,bins=bins,label='normal')
plt.hist(m2,bins=bins,label='abs(cauchy) <10',alpha=0.5)
plt.legend()
plt.show()
Vous pouvez voir que la grosse queue est forte. Considérez cela comme du bruit.
Ensuite, formez la somme des nombres aléatoires de la distribution de Koshi, et regardez la distribution de la dernière valeur.
from scipy.stats import cauchy
m1=[]
trial=1000000
for i in range(trial):
s = np.cumsum(np.random.normal(0, 1, 250))
m1.append(s[-1])
m2=[]
for i in range(trial):
s=[]
i=1
while i<=250:
s0 = (cauchy.rvs(0,1,1))
if abs(s0)<10:
s.append(s0)
i+=1
m2.append(np.cumsum(s)[-1])
bins=np.linspace(np.min(m2),np.max(m2),1000)
plt.hist(m1,bins=bins,label='normal')
plt.hist(m2,bins=bins,label='cauchy',alpha=0.5)
plt.legend()
plt.show()
Après tout, vous pouvez voir que l'écart est plus grand que la distribution normale.
Ensuite, créons une fonction AR1 qui utilise la distribution de Cauchy.
def ar1_c(rho,x0,c1,d,sigma,nsample):
e = cauchy.rvs(loc=0,scale=1,size=nsample)
c=c1-rho # constant of autoregressive mdoel
y=[]
y.append(e[0]+x0*rho+c+d)
for i,ee in enumerate(e[1:]):
y.append(rho*y[i]+ee+c+d*(i+1))
return y
Le test DF est effectué à l'aide de la série temporelle d'auto-retour du premier ordre créée.
def adfcheck_c(rho,x0,c1,d,sigma,nsample,trial):
prw=0
prwc=0
prwct=0
for i in range(trial):
y=ar1_c(rho,x0,c1,d,sigma,nsample)
rw=sm.tsa.adfuller(y,maxlag=0,regression='nc')[1]
rwc=sm.tsa.adfuller(y,maxlag=0,regression='c')[1]
rwct=sm.tsa.adfuller(y,maxlag=0,regression='ct')[1]
if rw>0.05:
prw+=1
if rwc>0.05:
prwc+=1
if rwct>0.05:
prwct+=1
print('nsample',nsample,prw/trial,prwc/trial,prwct/trial)
rho=0.90
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck_c(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.973 0.8775 0.8377
# nsample 20 0.9873 0.9588 0.9143
# nsample 250 0.8816 0.8425 0.0897
En raison de l'influence de la grosse queue, il s'agit dans ce cas d'une tendance temporelle avec dérive, mais lorsque la taille de l'échantillon est de 250, l'hypothèse nulle ne peut être rejetée au niveau de signification de 5%. Lorsque rho = 0,85, le résultat est le suivant.
nsample 7 0.9717 0.8801 0.8499 nsample 20 0.987 0.9507 0.903 nsample 250 0.8203 0.7057 0.0079
Lorsque la taille de l'échantillon est de 250, l'hypothèse nulle peut être rejetée au niveau de signification de 5%.
Ensuite, coupons la partie grosse queue de la distribution de Cauchy.
def cau(cut,nsample):
s=[]
i=1
while i<=nsample:
s0 = cauchy.rvs(0,1,1)
if abs(s0)<cut:
s.append(s0[0])
i+=1
return s
def ar1_c2(rho,x0,c1,d,sigma,nsample):
e = cau(10,nsample)
c=c1-rho # constant of autoregressive mdoel
y=[]
y.append(e[0]*sigma+x0*rho+c+d)
for i,ee in enumerate(e[1:]):
y.append(rho*y[i]+ee*sigma+c+d*(i+1))
return y
def adfcheck_c2(rho,x0,c1,d,sigma,nsample,trial):
prw=0
prwc=0
prwct=0
for i in range(trial):
y=ar1_c2(rho,x0,c1,d,sigma,nsample)
rw=sm.tsa.adfuller(y,maxlag=0,regression='nc')[1]
rwc=sm.tsa.adfuller(y,maxlag=0,regression='c')[1]
rwct=sm.tsa.adfuller(y,maxlag=0,regression='ct')[1]
if rw>0.05:
prw+=1
if rwc>0.05:
prwc+=1
if rwct>0.05:
prwct+=1
print('nsample',nsample,prw/trial,prwc/trial,prwct/trial)
rho=0.90
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck_c2(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.9969 0.9005 0.8462
# nsample 20 1.0 0.9912 0.9182
# nsample 250 1.0 1.0 0.0931
Le résultat était différent du précédent. Seules les tendances temporelles avec dérive peuvent rejeter l'hypothèse nulle, mais le niveau de signification est de 10%. Lorsque rho = 0,85
nsample 7 0.9967 0.9012 0.8466 nsample 20 1.0 0.9833 0.9121 nsample 250 1.0 1.0 0.0028
Le résultat est.
Ensuite, essayez le test Dicky Fuller étendu.
def adfcheck_c3(rho,x0,c1,d,sigma,nsample,trial):
prw=0
prwc=0
prwct=0
for i in range(trial):
y=ar1_c2(rho,x0,c1,d,sigma,nsample)
rw=sm.tsa.adfuller(y,regression='nc')[1]
rwc=sm.tsa.adfuller(y,regression='c')[1]
rwct=sm.tsa.adfuller(y,regression='ct')[1]
if rw>0.05:
prw+=1
if rwc>0.05:
prwc+=1
if rwct>0.05:
prwct+=1
print('nsample',nsample,prw/trial,prwc/trial,prwct/trial)
rho=0.90
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck_c3(rho,x0,c1,d,sigma,i,trial)
# nsample 7 0.94 0.8187 0.8469
# nsample 20 1.0 0.873 0.7969
# nsample 250 1.0 1.0 0.1253
Il ne peut pas être jugé correctement avec rho = 0,90. Si rho = 0,85, le résultat sera
nsample 7 0.9321 0.8109 0.8429 nsample 20 1.0 0.8512 0.7942 nsample 250 1.0 1.0 0.0252
J'ai pu rejeter la tendance dérive + temps.
Détection de la dérive et de la tendance temporelle lorsque l'écart type est nul ou faible
#dérive uniquement, l'écart type est nul
rho=0
c1=1.1
c1=(c1-1)/250+1
d=0
d=d/250
sigma=0#.4/np.sqrt(250)
trial=1#0000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
nsample 7 0.0 0.0 0.0 nsample 20 0.0 0.0 0.0 nsample 250 0.0 0.0 0.0
L'hypothèse nulle est rejetée du tout lorsqu'il n'y a qu'une dérive.
#dérive et tendance temporelle, l'écart type est nul
rho=0
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0#.4/np.sqrt(250)
trial=1#0000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
nsample 7 1.0 1.0 0.0 nsample 20 1.0 1.0 0.0 nsample 250 1.0 1.0 0.0 AR (1) + dérive + tendance temporelle rejette l'hypothèse nulle quelle que soit la taille de l'échantillon
Ensuite, ajoutez un écart type à la dérive.
# drift+ very small sigma
rho=0
c1=1.1
c1=(c1-1)/250+1
d=0
d=d/250
sigma=0.004#.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
nsample 7 0.9986 0.6379 0.7223 nsample 20 1.0 0.0396 0.1234 nsample 250 1.0 0.0 0.0
Ensuite, ajoutez un écart type à la tendance dérive + temps.
# drift+ very small sigma
rho=0
c1=1.1
c1=(c1-1)/250+1
d=1
d=d/250
sigma=0.004#.4/np.sqrt(250)
trial=10000
for i in [7,20,250]:
adfcheck(rho,x0,c1,d,sigma,i,trial)
nsample 7 1.0 0.9947 0.6877 nsample 20 1.0 1.0 0.1095 nsample 250 1.0 1.0 0.0
référence
statsmodels.tsa.stattools.adfuller
[Test Dicky-Fuller](https://ja.wikipedia.org/wiki/%E3%83%87%E3%82%A3%E3%83%83%E3%82%AD%E3%83%BC% E2% 80% 93% E3% 83% 95% E3% 83% A9% E3% 83% BC% E6% A4% 9C% E5% AE% 9A)
http://web.econ.ku.dk/metrics/Econometrics2_05_II/Slides/08_unitroottests_2pp.pdf