trajectory.jl 6.94 KB
Newer Older
1
2
3
4
5
6

abstract type AbstractTrajectory end
ContinuousObservations = AbstractVector{AbstractTrajectory}

struct Trajectory <: AbstractTrajectory
    m::ContinuousTimeModel
7
8
9
    values::Matrix{Int}
    times::Vector{Float64}
    transitions::Vector{Union{String,Nothing}}
10
11
end

12
# Operations
13
14
function +(σ1::AbstractTrajectory,σ2::AbstractTrajectory) end
function -(σ1::AbstractTrajectory,σ2::AbstractTrajectory) end
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140

# Top-level Lp distance function
"""
    `list_dist_lp(l_σ1, l_σ2; verbose, p, str_stat_list, str_stat_trajectory)`   

Function that computes Lp distance between two set of any dimensional trajectories.
...
# Arguments
- `l_σ1::Vector{AbstractTrajectory}` is the first set of trajectories. l_σ2 is the second.
- `verbose::Bool` If `true`, launch a verbose execution of the computation. 
- `str_stat_list::String` allows to set the statistic to apply on the distances of the trajectories. 
It is salso called linkage function in clustering field. Only "mean" is available for now on.
...
"""
function list_dist_lp(l_σ1::Vector{AbstractTrajectory}, l_σ2::Vector{AbstractTrajectory};  
                 verbose::Bool = false, p::Int = 1, str_stat_list::String = "mean", str_stat_trajectory::String = "mean")
    if str_stat_list == "mean"
        return dist_lp_mean(l_σ1, l_σ2; verbose = verbose, p = p, str_stat_trajectory = str_stat_trajectory)
    end
    error("Unrecognized statistic in dist_lp")
end
"""
    `dist_lp(σ1, σ2; verbose, p, str_stat)`   

Function that computes Lp distance between two trajectories of any dimension.
It computes Lp distances for each observed variable (contained in `get_obs_var(σ)`). and then apply a statistic on these distances.
Requires `get_obs_var(σ1) == get_obs_var(σ2)`, it is verified if they are simulated from the same model.
...
# Arguments
- `σ1::AbstractTrajectory` is the first trajectory. σ2 is the second.
- `verbose::Bool` If `true`, launch a verbose execution of the computation. 
- `str_stat::String` allows to set the statistic to apply on the distances computed  of the trajectories. 
Only "mean" is available for now on.
...
"""
function dist_lp(σ1::AbstractTrajectory, σ2::AbstractTrajectory;  
                 verbose::Bool = false, p::Int = 1, str_stat::String = "mean")
    if str_stat == "mean"
        return dist_lp_mean(σ1, σ2; verbose = verbose, p = p)
    end
    error("Unrecognized statistic in dist_lp")
end
"""
    `dist_lp(σ1, σ2, var; verbose, p, str_stat)`   

Function that computes Lp distance between two trajectories of any dimension.
It computes Lp distances for each observed variable (contained in `get_obs_var(σ)`). and then apply a statistic on these distances.
Requires `get_obs_var(σ1) == get_obs_var(σ2)`, it is verified if they are simulated from the same model.
...
# Arguments
- `σ1::AbstractTrajectory` is the first trajectory. σ2 is the second.
- `var::String` is an observed variable. Have to be contained in `get_obs_var(σ1)` and `get_obs_var(σ2)`.
- `verbose::Bool` If `true`, launch a verbose execution of the computation. 
...
"""
function dist_lp(σ1::AbstractTrajectory, σ2::AbstractTrajectory, var::String; 
                 verbose::Bool = false, p::Int = 1)
    if !is_bounded(σ1) || !is_bounded(σ2)
        @warn "Lp distance computed on unbounded trajectories. Result should be wrong"
    end
    return dist_lp(σ1[var], σ1["times"], σ2[var], σ2["times"])
end

# Distance function. Vectorized version
function dist_lp(x_obs::AbstractVector{Int}, t_x::Vector{Float64},  y_obs::AbstractVector{Int}, t_y::Vector{Float64}; 
                 verbose::Bool = false, p::Int = 1)
    current_y_obs = y_obs[1]
    current_t_y = t_y[2]
    idx = 1
    res = 0.0
    for i = 1:(length(x_obs)-1)
        if verbose
            @show i
            @show (t_x[i], x_obs[i])
            @show current_t_y
            @show t_x[i+1]
        end
        last_t_y = t_x[i]
        while current_t_y < t_x[i+1]
            rect =  abs(current_y_obs - x_obs[i]) * (current_t_y - last_t_y)
            res += rect
            if verbose
                println("-- in loop :")
                println("-- add : $rect abs($current_y_obs - $(x_obs[i])) * ($current_t_y - $(last_t_y)) / $(abs(current_y_obs - x_obs[i])) * $(current_t_y - last_t_y)")
                print("-- ")
                @show current_y_obs, current_t_y
            end
            idx += 1
            last_t_y = t_y[idx]
            current_y_obs = y_obs[idx]
            current_t_y = t_y[idx+1]
        end
        last_t_read = max(last_t_y, t_x[i])
        rect = abs(current_y_obs - x_obs[i]) * (t_x[i+1] - last_t_read)
        res += rect
        if verbose
            println("add : $rect abs($current_y_obs - $(x_obs[i])) * ($(t_x[i+1]) - $last_t_read) / $(abs(current_y_obs - x_obs[i])) * $(t_x[i+1] - last_t_read)")
            @show t_x[i+1], t_y[idx+1]
            @show y_obs[idx], x_obs[i]
            @show last_t_read
            @show current_y_obs
        end
    end
    return res
end
# For all the observed variables
function dist_lp_mean(σ1::AbstractTrajectory, σ2::AbstractTrajectory;  
                      verbose::Bool = false, p::Int = 1)
    if get_obs_var(σ1) != get_obs_var(σ2) error("Lp distances should be computed with the same observed variables") end
    res = 0.0
    for var in get_obs_var(σ1)
        res += dist_lp(σ1, σ2, var)
    end
    return res / length_obs_var(σ1)
end
# For a list of trajectories
function list_dist_lp_mean(l_σ1::Vector{AbstractTrajectory}, l_σ2::Vector{AbstractTrajectory};  
                 verbose::Bool = false, p::Int = 1, str_stat_trajectory::String = "mean")
    res = 0.0
    for σ1 in l_σ1
        for σ2 in l_σ2
            res += dist_lp(σ1, σ2; verbose = verbose, p = p, str_stat = str_stat_trajectory)
        end
    end
    return res / (length(l_σ1)*length(l_σ2))
end
141

142
# Properties of the trajectory
143
144
145
146
length_states(σ::AbstractTrajectory) = length(σ.times)
length_obs_var(σ::AbstractTrajectory) = size(σ.values)[2]
get_obs_var(σ::AbstractTrajectory) = (σ.m).g
is_bounded(σ::AbstractTrajectory) = σ.transitions[end] == nothing 
Bentriou Mahmoud's avatar
Bentriou Mahmoud committed
147

148
149
150
151
152
153
# Access to trajectory values
get_var_values(σ::AbstractTrajectory, var::String) = 
@view σ.values[:,(σ.m)._map_obs_var_idx[var]] 
get_state(σ::AbstractTrajectory, idx::Int) = @view σ.values[idx,:]
get_value(σ::AbstractTrajectory, var::String, idx::Int) = 
σ.values[idx,(σ.m)._map_obs_var_idx[var]] 
154
function δ(σ1::AbstractTrajectory,t::Float64) end
155

156
157
158
# Get var values ["I"]
function getindex(σ::AbstractTrajectory, var::String)
    if var  == "times"
159
        return σ.times
160
    elseif var == "transitions"
161
162
        return σ.transitions
    else
163
164
165
166
        return get_var_values(σ, var)
    end
end
# Get i-th state [i]
167
getindex(σ::AbstractTrajectory, idx::Int) = get_state(σ, idx)
168
# Get i-th value of var ["I", idx]
169
function getindex(σ::AbstractTrajectory, var::String, idx::Int)
170
171
172
173
174
175
    if var  == "times"
        return σ.times[idx]
    elseif var == "transitions"
        return σ.transitions[idx]
    else
        return get_value(σ, var, idx)
176
177
    end
end
178