Accessing simulation data

A single run generates a SamplePath object which is a recursive array provided by RecursiveArrayTools.jl.

sample_path = simulate(state, model, SortingDirect(), tfinal = 100.0, save_points = 0:25:100.0)

sample_path
t: 5-element Array{Float64,1}:
   0.0
  25.0
  50.0
  75.0
 100.0
x: 5-element Array{Array{Int64,1},1}:
 [10, 0, 0, 0, 0]  
 [3, 7, 1, 9, 61]  
 [1, 9, 0, 13, 102]
 [1, 9, 0, 11, 101]
 [1, 9, 1, 12, 118]

Access the state at the first time point:

sample_path[1]
5-element Array{Int64,1}:
 10
  0
  0
  0
  0

Access the third component at the second time point:

sample_path[3,2]
1

Access entire history for the third and fourth components:

sample_path[[3,4], :]
2×5 Array{Int64,2}:
 0  1   0   0   1
 0  9  13  11  12

Running multiple simulations generates an Ensemble, a collection of SamplePath objects; that is, Vector{SamplePath}. These objects also support indexing:

# ensemble of 10 trajectories
ensemble = [simulate(state, model, SortingDirect(), tfinal = 100.0, save_points = 0:25:100.0) for _ in 1:10];
10-element Array{SamplePath{Int64,2,Array{Array{Int64,1},1},Array{Float64,1}},1}:
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [4, 6, 2, 12, 96], [0, 10, 1, 36, 265], [1, 9, 0, 30, 277], [0, 10, 0, 25, 287]]  
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [0, 10, 1, 17, 159], [0, 10, 0, 15, 159], [2, 8, 0, 13, 164], [0, 10, 0, 15, 158]]
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [9, 1, 0, 1, 6], [1, 9, 1, 13, 114], [1, 9, 0, 11, 114], [0, 10, 1, 20, 187]]     
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [7, 3, 0, 1, 31], [2, 8, 0, 3, 24], [1, 9, 0, 12, 49], [3, 7, 0, 11, 51]]         
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [3, 7, 0, 2, 9], [7, 3, 0, 5, 11], [6, 4, 0, 4, 10], [7, 3, 2, 12, 22]]           
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [2, 8, 1, 10, 47], [0, 10, 0, 15, 116], [1, 9, 0, 24, 111], [0, 10, 0, 11, 114]]  
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [6, 4, 0, 3, 1], [3, 7, 0, 9, 33], [2, 8, 0, 10, 31], [5, 5, 0, 11, 33]]          
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [0, 10, 1, 13, 55], [1, 9, 1, 21, 181], [2, 8, 0, 20, 248], [0, 10, 0, 23, 241]]  
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [0, 10, 1, 13, 158], [1, 9, 0, 16, 160], [0, 10, 0, 11, 160], [1, 9, 0, 9, 160]]  
 t: [0.0, 25.0, 50.0, 75.0, 100.0]
x: Array{Int64,1}[[10, 0, 0, 0, 0], [0, 10, 1, 14, 102], [0, 10, 0, 18, 133], [1, 9, 0, 13, 135], [0, 10, 0, 10, 132]]

Retrieve an individual sample path:

ensemble[2]
t: 5-element Array{Float64,1}:
   0.0
  25.0
  50.0
  75.0
 100.0
x: 5-element Array{Array{Int64,1},1}:
 [10, 0, 0, 0, 0]   
 [0, 10, 1, 17, 159]
 [0, 10, 0, 15, 159]
 [2, 8, 0, 13, 164] 
 [0, 10, 0, 15, 158]

Index into a SamplePath:

ensemble[2][[3,4], :]
2×5 Array{Int64,2}:
 0   1   0   0   0
 0  17  15  13  15

Working with DataFrames

Both SamplePath and Ensemble implement the IterableTables.jl interface. This means that they act as data sources that can be converted into any supported data sink. For example, you can convert an Ensemble into a DataFrame:

using DataFrames

DataFrame(ensemble)

50 rows × 7 columns

trialtX1X2X3X4X5
Int64Float64Int64Int64Int64Int64Int64
110.0100000
2125.04621296
3150.0010136265
4175.019030277
51100.0010025287
620.0100000
7225.0010117159
8250.0010015159
9275.028013164
102100.0010015158
1130.0100000
12325.091016
13350.019113114
14375.019011114
153100.0010120187
1640.0100000
17425.0730131
18450.0280324
19475.01901249
204100.03701151
2150.0100000
22525.037029
23550.0730511
24575.0640410
255100.07321222
2660.0100000
27625.02811047
28650.0010015116
29675.019024111
306100.0010011114
3170.0100000
32725.064031
33750.0370933
34775.02801031
357100.05501133
3680.0100000
37825.001011355
38850.019121181
39875.028020248
408100.0010023241
4190.0100000
42925.0010113158
43950.019016160
44975.0010011160
459100.01909160
46100.0100000
471025.0010114102
481050.0010018133
491075.019013135
5010100.0010010132

Each component of the state vector appears as a column, along with trajectory trial and timestamp t columns. If the conversion does not produce the expected result, one may be able to force the correct behavior using the tablefy function:

import BioSimulator: tablefy

DataFrame(tablefy(sample_path)) # DataFrame(sample_path) is currently incorrect

5 rows × 6 columns

tX1X2X3X4X5
Float64Int64Int64Int64Int64Int64
10.0100000
225.0371961
350.019013102
475.019011101
5100.019112118

Saving simulation results

Because SamplePath and Ensemble support iteration, you can save simulation data directly using Julia's I/O interface.

Example here

The easiest approach takes advantage of IterableTables.jl:

using CSV

CSV.write(file, result, delim = '\t')

Obtaining summary statistics

Summary statistics for ensembles are supported directly:

using Statistics

mean(ensemble)
std(result)
var(result)