Tau-hoppande

Inom sannolikhetsteorin är tau-leaping , eller τ-leaping , en ungefärlig metod för simulering av ett stokastiskt system . Den är baserad på Gillespie-algoritmen , som utför alla reaktioner under ett intervall av längd tau innan uppdatering av benägenhetsfunktionerna. Genom att uppdatera priserna mindre ofta möjliggör detta ibland effektivare simulering och därmed hänsyn till större system.

Många varianter av den grundläggande algoritmen har övervägts.

Algoritm

Algoritmen är analog med Eulermetoden för deterministiska system, men istället för att göra en fast förändring

$x(t+\tau )=x(t)+\tau x'(t)$

förändringen är

$x(t+\tau )=x(t)+P(\tau x'(t))$

där $P(\tau x'(t))$ är en Poisson- fördelad slumpvariabel med medelvärdet $\tau x'(t)$ .

Givet ett tillstånd $\mathbf {x} (t)=\{X_{i}(t)\}$ med händelser $E_{j}$ förekommer med hastigheten $R_{j}(\mathbf {x} (t))$ och med tillståndsändringsvektorer $\mathbf {v} _{ij}$ (där $i$ indexerar tillståndsvariablerna och $j$ indexerar händelserna), är metoden som följer:

1. Initiera modellen med initiala villkor $\mathbf {x} (t_{0})=\{X_{i}(t_{0})\}$ .
2. Beräkna händelsehastigheterna $R_{j}(\mathbf {x} (t))$ .
3. Välj ett tidssteg $\tau$ . Detta kan vara fixat eller med någon algoritm beroende på de olika händelsehastigheterna.
4. För varje händelse genererar $E_{j}$${\displaystyle K_{j}\sim {\text{Poisson}}(R_{j}\tau )} ,$ som är antalet gånger varje händelse inträffar under tidsintervallet $[t,t+\tau )$ .
5. Uppdatera tillståndet med

   $\mathbf {x} (t+\tau )=\mathbf {x} (t)+\summa _{ j}K_{j}v_{ij}$

   där $v_{ij}$ är förändringen av tillståndsvariabel $X_{i}$ på grund av händelse $E_{j }$ . Vid denna tidpunkt kan det vara nödvändigt att kontrollera att inga populationer har nått orealistiska värden (som att en population blir negativ på grund av Poisson-variabeln K j {\displaystyle K_{j $)$ .

6. Upprepa från steg 2 och framåt tills något önskat villkor är uppfyllt (t.ex. en viss tillståndsvariabel når 0, eller tiden $t_{1}$ uppnås).

Algoritm för effektivt val av stegstorlek

Denna algoritm beskrivs av Cao et al. Tanken är att begränsa den relativa förändringen i varje händelsehastighet $R_{j}$ med en specificerad tolerans $\epsilon$ (Cao et al. rekommenderar $\epsilon =0,03$ , även om det kan bero på modellspecifikationer). Detta uppnås genom att begränsa den relativa förändringen i varje tillståndsvariabel $X_{i}$ med ${\displaystyle \epsilon /g_{i}} ,$ där $g_{i}$ beror på den hastighet som ändras mest för en given förändring i $X_{i}$ . Typiskt $g_{i}$ lika med den högsta ordningens händelsefrekvens, men detta kan vara mer komplext i olika situationer (särskilt epidemiologiska modeller med icke-linjära händelsefrekvenser).

Denna algoritm kräver vanligtvis beräkning av $2N$ hjälpvärden (där $N$ är antalet tillståndsvariabler $X_{i}$ ), och bör endast kräva återanvändning av tidigare beräknade värden $R_{j}(\mathbf {x} )$ . En viktig faktor i detta eftersom $X_{i}$ är ett heltalsvärde, då finns det ett minimivärde med vilket det kan ändras, vilket förhindrar att den relativa förändringen i $R_{j}$ begränsas med 0, vilket skulle resultera i att $\tau$ också tenderar till 0.

1. För varje tillståndsvariabel $X_{i}$ beräknar du hjälpvärdena

   $\mu _{i}(\mathbf {x} )=\summa _{j}v_{ij}R_{j}(\mathbf {x} )$

   $\sigma _{i }^{2}(\mathbf {x} )=\summa _{j}v_{ij}^{2}R_{j}(\mathbf {x} )$

2. För varje tillståndsvariabel $X_{i}$ , bestäm den högsta ordningens händelse som den är involverad i och erhåll $g_{i}$
3. Beräkna tidssteg $\tau$ som

   $\tau =\min _{i}{\left\{{\frac {\max {\{\epsilon X_{ i}/g_{i},1\}}}{|\mu _{i}(\mathbf {x} )|}},{\frac {\max {\{\epsilon X_{i}/g_{ i},1\}}^{2}}{\sigma _{i}^{2}(\mathbf {x} )}}\right\}}$

Denna beräknade $\tau$ används sedan i steg 3 i $\tau$ hoppalgoritmen.