Merge pull request #275 from Neuroblox/scopetype

merging now with the understanding that we need to go back to the issue of parameter sharing.  blox_utilities.jl is all over the place and needs to be sorted out.  Good news is that Spectral DCM is working with AD.

Project.toml

@@ -6,6 +6,7 @@ version = "0.1.0"
 [deps]
 AbstractFFTs = "621f4979-c628-5d54-868e-fcf4e3e8185c"
 Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
+ChainRules = "082447d4-558c-5d27-93f4-14fc19e9eca2"
 Combinatorics = "861a8166-3701-5b0c-9a16-15d98fcdc6aa"
 DSP = "717857b8-e6f2-59f4-9121-6e50c889abd2"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
@@ -40,6 +41,7 @@ Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+RuntimeGeneratedFunctions = "7e49a35a-f44a-4d26-94aa-eba1b4ca6b47"
 SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
 Serialization = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
 SignalAnalysis = "df1fea92-c066-49dd-8b36-eace3378ea47"
@@ -119,4 +121,4 @@ Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c"
 Turing = "fce5fe82-541a-59a6-adf8-730c64b5f9a0"
 
 [targets]
-test = ["DifferentialEquations","LinearAlgebra", "ForwardDiff", "Distributions", "Flux", "Random", "SparseArrays", "Optim", "Optimization", "OptimizationOptimJL", "OptimizationOptimisers", "OrdinaryDiffEq", "StochasticDiffEq", "SafeTestsets", "Turing", "Test", "Tracker"]
+test = ["DifferentialEquations", "LinearAlgebra", "ForwardDiff", "Distributions", "Flux", "Random", "SparseArrays", "Optim", "Optimization", "OptimizationOptimJL", "OptimizationOptimisers", "OrdinaryDiffEq", "StochasticDiffEq", "SafeTestsets", "Turing", "Test", "Tracker"]

src/Neuroblox.jl

@@ -55,13 +55,16 @@ abstract type BloxConnectComplex <: BloxConnection end
 abstract type BloxConnectMultiFloat <: BloxConnection end
 abstract type BloxConnectMultiComplex <: BloxConnection end
 
+# dictionary type for Blox parameters
+Para_dict = Dict{Symbol, Union{<: Real, Num}}
+
 include("Neurographs.jl")
+include("blox/blox_utilities.jl")
 include("utilities/spectral_tools.jl")
 include("utilities/learning_tools.jl")
 include("control/controlerror.jl")
 include("measurementmodels/fmri.jl")
 include("datafitting/spectralDCM.jl")
-include("blox/blox_utilities.jl")
 include("blox/neural_mass.jl")
 include("blox/cortical_blox.jl")
 include("blox/canonicalmicrocircuit.jl")

src/blox/blox_utilities.jl

@@ -1,11 +1,58 @@
-function scope_dict(para_dict::Dict{Symbol, Union{Real,Num}})
-    para_dict_copy = copy(para_dict)
-    for (n,v) in para_dict_copy
-        if typeof(v) == Num
-            para_dict_copy[n] = ParentScope(v)
+function progress_scope(params; lvl=0)
+    para_list = []
+    for p in params
+        pp = ModelingToolkit.unwrap(p)
+        if ModelingToolkit.hasdefault(pp)
+            d = ModelingToolkit.getdefault(pp)
+            if typeof(d)==SymbolicUtils.BasicSymbolic{Real}
+                if lvl==0
+                    pp = ParentScope(pp)
+                else
+                    pp = DelayParentScope(pp,lvl)
+                end
+            end
+        end
+        push!(para_list,ModelingToolkit.wrap(pp))
+    end
+    return para_list
+end
+
+"""
+This function progresses the scope of parameters and leaves floating point values untouched
+"""
+function progress_scope(args...)
+    paramlist = []
+    for p in args
+        if p isa Float64
+            push!(paramlist, p)
         else
-            para_dict_copy[n] = (@parameters $n=v)[1]
+            p = ParentScope(p)
+            # pp = ModelingToolkit.unwrap(p)
+            # if ModelingToolkit.hasdefault(pp)
+            #     d = ModelingToolkit.getdefault(pp)
+            #     if typeof(d)==SymbolicUtils.BasicSymbolic{Real}
+            #         pp = ParentScope(pp)
+            #     end
+            # end
+            # push!(para_list,ModelingToolkit.wrap(pp))
+            push!(paramlist, p)
         end
     end
-    return para_dict_copy
+    return paramlist
 end
+
+"""
+    This function compiles already existing parameters with floats after making them parameters.
+    Keyword arguments are used because parameter definition requires names, not just values
+"""
+function compileparameterlist(;kwargs...)
+    paramlist = []
+    for (kw, v) in kwargs
+        if v isa Float64
+            paramlist = vcat(paramlist, @parameters $kw = v)
+        else
+            paramlist = vcat(paramlist, v)
+        end
+    end
+    return paramlist
+end

src/blox/neural_mass.jl

@@ -14,23 +14,19 @@ end
 
 mutable struct HarmonicOscillatorBlox <: NeuralMassBlox
     # all parameters are Num as to allow symbolic expressions
-    p_dict::Dict{Symbol,Union{Real,Num}}
     connector::Num
     noDetail::Vector{Num}
     detail::Vector{Num}
     initial::Dict{Num, Tuple{Float64, Float64}}
     odesystem::ODESystem
     function HarmonicOscillatorBlox(;name, ω=25*(2*pi), ζ=1.0, k=625*(2*pi), h=35.0)
-        para_dict = scope_dict(Dict{Symbol,Union{Real,Num}}(:ω => ω,:ζ => ζ,:k => k,:h => h))
-        ω=para_dict[:ω]
-        ζ=para_dict[:ζ]
-        k=para_dict[:k]
-        h=para_dict[:h]
+        params = progress_scope(@parameters ω=ω ζ=ζ k=k h=h)
         sts    = @variables x(t)=1.0 y(t)=1.0 jcn(t)=0.0
+        ω, ζ, k, h = params
         eqs    = [D(x) ~ y-(2*ω*ζ*x)+ k*(2/π)*(atan((jcn)/h))
                   D(y) ~ -(ω^2)*x]
-        odesys = ODESystem(eqs, t, sts, values(para_dict); name=name)
-        new(para_dict, odesys.x,[odesys.x],[odesys.x,odesys.y],
+        odesys = ODESystem(eqs, t, sts, params; name=name)
+        new(odesys.x,[odesys.x],[odesys.x,odesys.y],
             Dict(odesys.x => (-1.0,1.0), odesys.y => (-1.0,1.0)),
             odesys)
     end
@@ -43,23 +39,19 @@ const harmonic_oscillator = HarmonicOscillatorBlox
 #    return HarmonicOscillatorImage
 
 mutable struct JansenRitCBlox <: NeuralMassBlox
-    p_dict::Dict{Symbol,Union{Real,Num}}
     connector::Num
     noDetail::Vector{Num}
     detail::Vector{Num}
     initial::Dict{Num, Tuple{Float64, Float64}}
     odesystem::ODESystem
     function JansenRitCBlox(;name, τ=0.001, H=20.0, λ=5.0, r=0.15)
-        para_dict = scope_dict(Dict{Symbol,Union{Real,Num}}(:τ => τ,:H => H,:λ => λ,:r => r))
-        τ=para_dict[:τ]
-        H=para_dict[:H]
-        λ=para_dict[:λ]
-        r=para_dict[:r]
+        params = progress_scope(@parameters τ=τ H=H λ=λ r=r)
         sts    = @variables x(t)=1.0 y(t)=1.0 jcn(t)=0.0
+        τ, H, λ, r = params
         eqs    = [D(x) ~ y - ((2/τ)*x),
                 D(y) ~ -x/(τ*τ) + (H/τ)*((2*λ)/(1 + exp(-r*(jcn))) - λ)]
-        odesys = ODESystem(eqs, t, sts, values(para_dict); name=name)
-        new(para_dict, odesys.x,[odesys.x],[odesys.x,odesys.y],
+        odesys = ODESystem(eqs, t, sts, params; name=name)
+        new(odesys.x,[odesys.x],[odesys.x,odesys.y],
             Dict(odesys.x => (-1.0,1.0), odesys.y => (-1.0,1.0)),
             odesys)
     end
@@ -68,23 +60,19 @@ end
 const jansen_ritC = JansenRitCBlox
 
 mutable struct  JansenRitSCBlox <: NeuralMassBlox
-    p_dict::Dict{Symbol,Union{Real,Num}}
     connector::Num
     noDetail::Vector{Num}
     detail::Vector{Num}
     initial::Dict{Num, Tuple{Float64, Float64}}
     odesystem::ODESystem
     function JansenRitSCBlox(;name, τ=0.014, H=20.0, λ=400.0, r=0.1)
-        para_dict = scope_dict(Dict{Symbol,Union{Real,Num}}(:τ => τ,:H => H,:λ => λ,:r => r))
-        τ=para_dict[:τ]
-        H=para_dict[:H]
-        λ=para_dict[:λ]
-        r=para_dict[:r]
+        params = progress_scope(@parameters τ=τ H=H λ=λ r=r)
         sts    = @variables x(t)=1.0 y(t)=1.0 jcn(t)=0.0
+        τ, H, λ, r = params
         eqs    = [D(x) ~ y - ((2/τ)*x),
                   D(y) ~ -x/(τ*τ) + (H/τ)*((2*λ)/(1 + exp(-r*(jcn))) - λ)]
-        odesys = ODESystem(eqs, t, sts, values(para_dict); name=name)
-        new(para_dict, odesys.x,[odesys.x],[odesys.x,odesys.y],
+        odesys = ODESystem(eqs, t, sts, params; name=name)
+        new(odesys.x,[odesys.x],[odesys.x,odesys.y],
             Dict(odesys.x => (-1.0,1.0), odesys.y => (-1.0,1.0)),
             odesys)
     end