ObservableCollector.jl

Easily define observables to be measured.

Introduction

Simulating a physical/biological/chemical,... system typically involes the following steps

Setup a state with some initial conditions.
Evole state through some deterministic or stochastic dynamics for a number of steps.
Measure a host of observables; possibly only if some condition is fulfilled.
Go back to 2. until the simulation reaches a halting condition.

Step 3. typically involves preparing a vector for each quantity to be measured (or a dictionary) beforehand, and then pushing the measurements onto them.

The package at hand implements a convinience layer around step 3.

Terminology:

state The current state of the model. Can be any data type or user defined structure. Contains all necessary information to determine observables from.
observable Function of state and possibly simulation time. In the an observable will refer to a function f(state, time)::Any.

Let us look at a quick (and mostly trivial) example. The state is simply a fixed-size vector. In each timestep we add a bit of random noise to every entry (effectively simulating many random walks on the real line). We'd like to measure mean and variance over time. Additionally we imagine the steps happening at rate 1.0 and thus draw the waiting-time between steps from an exponential distribution.

using ObservableCollector
using Distributions

p = Gamma(100, 1.0) # sum of exponentially distr. rvs is gamma-distributed.
state = fill(0.0, 1000)

obs! = @observations begin
  @condition (s,t)->t%100==0 begin
    "step" --> (s,t)->t
    "dt" --> (s,t)->rand(p)
    "m" --> (s,t)->mean(s)
    "var" --> (s,t)->var(s)
  end
end # note the ! in the name I chose. `obs!` will mutate its first argument.

X = [] # empty array to store all measured quantities in.
for t in 0:10^5 # make 10^5 steps
  obs!(X, state, t)
  state .+= 0.05*randn(length(state))
end

X

4004-element Array{Any,1}:
 (name = :step, val = 0)                 
 (name = :dt, val = 83.97155298518244)   
 (name = :m, val = 0.0)                  
 (name = :var, val = 0.0)                
 (name = :step, val = 100)               
 (name = :dt, val = 100.00740230972015)  
 (name = :m, val = 0.0003930655995079922)
 (name = :var, val = 0.26160886476073925)
 (name = :step, val = 200)               
 (name = :dt, val = 96.85494971173976)   
 ⋮                                       
 (name = :var, val = 255.93843105199122) 
 (name = :step, val = 99900)             
 (name = :dt, val = 93.37071123687382)   
 (name = :m, val = -0.5043983754920668)  
 (name = :var, val = 256.20664878355996) 
 (name = :step, val = 100000)            
 (name = :dt, val = 99.27002411873116)   
 (name = :m, val = -0.4984856006539832)  
 (name = :var, val = 256.795375882665)

This storage format is not very appealing. It is simply a linear record of observables collected throughout the run. We may use the timeseries function to cast it into a more useful form

timeseries(X)

Dict{Symbol,Any} with 4 entries:
  :m    => Union{Missing, Float64}[0.0, 0.000393066, -0.00328986, 0.00407247, -…
  :step => Union{Missing, Int64}[0, 100, 200, 300, 400, 500, 600, 700, 800, 900…
  :var  => Union{Missing, Float64}[0.0, 0.261609, 0.496245, 0.717733, 0.942514,…
  :dt   => Union{Missing, Float64}[83.9716, 100.007, 96.8549, 102.142, 111.19, …

Hints

t needn't necessarily be time. You may use it as some auxilliary variable, or even as a boolean switch to turn collection on/off. It also need not be a scalar, but could e.g. be a tuple.
Observables are (by my definition) functions of state only and thus have no access to the measurement history. If you require observables to depend on previous measurements, store those in the state. (TODO: example)

Available methods

# ObservableCollector.@observations — Macro.

@observations(block)

Creates an anonymous function with arguments (output::Vector{Any}, state, time) that collects observables under conditions defined in block.

Example:

output = []
obs = @observations begin
    @condition (s,t)->t>=3 begin
        "A" --> (s,t)->s
        "T" --> (s,t)->t
    end
end
for t in 1:4
    obs(output,1,t)
end

julia> @show output
output = Any[(name = :A, val = 1), (name = :T, val = 3), (name = :A, val = 1), (name = :T, val = 4)]
4-element Array{Any,1}:
 (name = :A, val = 1)
 (name = :T, val = 3)
 (name = :A, val = 1)
 (name = :T, val = 4)

source

# ObservableCollector.@condition — Macro.

@condition(cond, block)

Define multiple observables to be evaluated under the same condition. Use --> to define pairs of name --> map.

Example

@condition (s,t)->t>=3 begin
    "S" --> (s,t)->s^2
    "T" --> (s,t)->t
end

source

@condition(cond, name, map)

Creates a clause to collect an observable map(state,time) called name under condition cond(state,time).