In RETURN-project/BenchmarkRecovery: Benchmark recovery indicators derived from remote sensing time series

# Default chunk configuration
knitr::opts_chunk$set(echo = TRUE, 
                      eval = TRUE, 
                      fig.width=8, 
                      fig.height=6) 

# Load libraries
library(BenchmarkRecovery)

Comparing decay functions

# Set the problem
## Set the times
ts <- seq(10, 100, by = 1)

## Set the parameters (the same for all)
offset <- 0
pert <- 5
tpert <- 15
thalf <- 3
noise <- 0.1

# Generate the time series
ys_p <- piecewise(ts, offset, pert, tpert, thalf, noise)
ys_e <- exponential(ts, offset, pert, tpert, thalf, noise)
ys_r <- realistic(ts, offset, pert, tpert, thalf, noise)

Plots

Complete time series

# Plot the time series
plot(ts, ys_p, 'l', col = 'black', xlab = 'time', ylab = 'state')
lines(ts, ys_e, col = 'red')
lines(ts, ys_r, col = 'blue')
legend('topright', legend = c('piecewise', 'exponential', 'realistic'),
       col = c('black', 'red', 'blue'), lty = 1, cex = 1)

Detail (perturbation)

# Plot the time series
plot(ts, ys_p, 'l', col = 'black', xlab = 'time', ylab = 'state', xlim = c(tpert - thalf, tpert + thalf))
lines(ts, ys_e, col = 'red')
lines(ts, ys_r, col = 'blue')
legend('topright', legend = c('piecewise', 'exponential', 'realistic'),
       col = c('black', 'red', 'blue'), lty = 1, cex = 1)

Detail (transient)

# Plot the time series
plot(ts, ys_p, 'l', col = 'black', xlab = 'time', ylab = 'state', xlim = c(tpert, tpert + 7*thalf))
lines(ts, ys_e, col = 'red')
lines(ts, ys_r, col = 'blue')
legend('topright', legend = c('piecewise', 'exponential', 'realistic'),
       col = c('black', 'red', 'blue'), lty = 1, cex = 1)

Detail (asymptotics / long term dynamics)

# Plot the time series
plot(ts, ys_p, 'l', col = 'black', xlab = 'time', ylab = 'state', xlim = c(0.5*max(ts), max(ts)), ylim = c(-4*noise, 4*noise))
lines(ts, ys_e, col = 'red')
lines(ts, ys_r, col = 'blue')
legend('topright', legend = c('piecewise', 'exponential', 'realistic'),
       col = c('black', 'red', 'blue'), lty = 1, cex = 1)

Compare standard deviations

sds <- c(sd(ys_p), sd(ys_e), sd(ys_r))
print(sds)

RETURN-project/BenchmarkRecovery documentation built on July 13, 2021, 5:47 p.m.