Langevin1D {Langevin} | R Documentation |
Langevin1D
calculates the Drift and Diffusion vectors (with errors)
for a given time series.
Langevin1D(data, bins, steps, sf = ifelse(is.ts(data), frequency(data), 1), bin_min = 100, reqThreads = -1)
data |
a vector containing the time series or a time-series object. |
bins |
a scalar denoting the number of |
steps |
a vector giving the τ steps to calculate the conditional moments (in samples (=τ * sf)). |
sf |
a scalar denoting the sampling frequency (optional if |
bin_min |
a scalar denoting the minimal number of events per |
reqThreads |
a scalar denoting how many threads to use. Defaults to
|
Langevin1D
returns a list with thirteen components:
D1 |
a vector of the Drift coefficient for each |
eD1 |
a vector of the error of the Drift coefficient for each
|
D2 |
a vector of the Diffusion coefficient for each |
eD2 |
a vector of the error of the Driffusion coefficient for
each |
D4 |
a vector of the fourth Kramers-Moyal coefficient for each
|
mean_bin |
a vector of the mean value per |
density |
a vector of the number of events per |
M1 |
a matrix of the first conditional moment for each
τ. Rows corespond to |
eM1 |
a matrix of the error of the first conditional moment
for each τ. Rows corespond to |
M2 |
a matrix of the second conditional moment for each
τ. Rows corespond to |
eM2 |
a matrix of the error of the second conditional moment
for each τ. Rows corespond to |
M4 |
a matrix of the fourth conditional moment for each
τ. Rows corespond to |
U |
a vector of the |
Philip Rinn
# Set number of bins, steps and the sampling frequency bins <- 20; steps <- c(1:5); sf <- 1000; #### Linear drift, constant diffusion # Generate a time series with linear D^1 = -x and constant D^2 = 1 x <- timeseries1D(N=1e6, d11=-1, d20=1, sf=sf); # Do the analysis est <- Langevin1D(x, bins, steps, sf, reqThreads=2); # Plot the result and add the theoretical expectation as red line plot(est$mean_bin, est$D1); lines(est$mean_bin, -est$mean_bin, col='red'); plot(est$mean_bin, est$D2); abline(h=1, col='red'); #### Cubic drift, constant diffusion # Generate a time series with cubic D^1 = x - x^3 and constant D^2 = 1 x <- timeseries1D(N=1e6, d13=-1, d11=1, d20=1, sf=sf); # Do the analysis est <- Langevin1D(x, bins, steps, sf, reqThreads=2); # Plot the result and add the theoretical expectation as red line plot(est$mean_bin, est$D1); lines(est$mean_bin, est$mean_bin - est$mean_bin^3, col='red'); plot(est$mean_bin, est$D2); abline(h=1, col='red');