Fehlerberechnung für eine Zeitreihe mittels AutokorrelationPython

Python-Programme
Anonymous
 Fehlerberechnung für eine Zeitreihe mittels Autokorrelation

Post by Anonymous »

Ich versuche, die Autokorrelationszeit für eine Zeitreihe zu berechnen, die aus molekulardynamischen Trajektorienrahmen generiert wurde. Ich habe mehrere Methoden zur Schätzung der Autokorrelationszeit verwendet, da ich sie zur Berechnung des statistischen Fehlers einer Observablen benötige. Mein erster Versuch war die Blockmittelung. Das Blockmittelungsdiagramm für den Standardfehler des Mittelwerts (SEM) zeigt kein klares Plateau. Der Code, den ich für die Blockmittelung verwendet habe:

Code: Select all

sem=[]

#return sems
for i in range(1,1000):
block_size.append(i)
block_mean=[]
no_blocks=0
for l in range(0, len(nums),i):
if len(nums[l:l+i])==i:
block_mean.append(np.mean(nums[l:l+i]))
no_blocks=no_blocks+1
sem.append(np.std(block_mean)/np.sqrt(no_blocks))

Beispielausgabe (Blockmittelung SEM vs. Blockgröße):
Image
Meine Frage: Wenn es im SEM aufgrund der Blockmittelung kein klares Plateau gibt, ist es dann trotzdem sinnvoll, den statistischen Fehler für diese Daten abzuschätzen? Wenn ja, welche Methoden/Vorsichtsmaßnahmen sollte ich verwenden, um eine zuverlässige Fehlerschätzung (oder eine konservative Grenze) aus einer einzelnen langen Trajektorie von Molecular Dynamics-Frames zu erhalten?
Ich habe versucht, andere Methoden wie statsmodel.tsa.stattools.acf für die Berechnung der Autokorrelationszeit zu verwenden. Der Fehler, den ich am Ende erhalte, ist 0,0008712, was meiner Meinung nach etwas unrealistisch erscheint. Gibt es eine andere Möglichkeit, einen genaueren Fehler für diese Daten zu berechnen?
Ich habe den folgenden Code verwendet, um den Fehler mithilfe von statsmodel.tsa.stattools.acf zu berechnen,

Code: Select all

def tau_cal(acf_out,nums):
positive_acf =np.where(acf_out > 0.0000)[0]

result = [acf_out[x] for x in positive_acf]
tau = 1 + 2 * np.sum(result[1:-1])
Neff = len(nums)/tau
error = (statistics.stdev(nums))/np.sqrt(Neff)
return error
time = statsmodel.tsa.stattools.acf(nums)
error = tau_cal(acf_out,nums)

Quick Reply

Change Text Case: 
   
  • Similar Topics
    Replies
    Views
    Last post