diff --git a/dev/.documenter-siteinfo.json b/dev/.documenter-siteinfo.json index 8a5aaa4..476cc6c 100644 --- a/dev/.documenter-siteinfo.json +++ b/dev/.documenter-siteinfo.json @@ -1 +1 @@ -{"documenter":{"julia_version":"1.10.4","generation_timestamp":"2024-07-24T21:38:20","documenter_version":"1.5.0"}} \ No newline at end of file +{"documenter":{"julia_version":"1.10.4","generation_timestamp":"2024-07-26T18:27:18","documenter_version":"1.5.0"}} \ No newline at end of file diff --git a/dev/api/index.html b/dev/api/index.html index 0ccd9c4..e289f6f 100644 --- a/dev/api/index.html +++ b/dev/api/index.html @@ -1,5 +1,5 @@ -API · Neuroblox
GitHub

API Documentation

Neuroblox.BalloonModel — Type

Arguments:

  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • lnκ: logarithmic prefactor to signal decay H[1], set to 0 for standard parameter value.
  • lnτ: logarithmic prefactor to transit time H[3], set to 0 for standard parameter value.
  • lnϵ: logarithm of ratio of intra- to extra-vascular signal

NB: the prefix ln of the variables u, ν, q as well as the parameters κ, τ denotes their transformation into logarithmic space to enforce their positivity. This transformation is considered in the derivates of the model equations below.

Citations:

  1. Stephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040
  2. Hofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407
source
Neuroblox.Generic2dOscillator — Type
Generic2dOscillator(name, namespace, ...)
+API · Neuroblox
GitHub
API Documentation
Neuroblox.BalloonModel — Type
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • lnκ: logarithmic prefactor to signal decay H[1], set to 0 for standard parameter value.
  • lnτ: logarithmic prefactor to transit time H[3], set to 0 for standard parameter value.
  • lnϵ: logarithm of ratio of intra- to extra-vascular signal
NB: the prefix ln of the variables u, ν, q as well as the parameters κ, τ denotes their transformation into logarithmic space to enforce their positivity. This transformation is considered in the derivates of the model equations below. 
Citations:
  1. Stephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040
  2. Hofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407
source
Neuroblox.Generic2dOscillator — Type
Generic2dOscillator(name, namespace, ...)
 
 The Generic2dOscillator model is a generic dynamic system with two state
 variables. The dynamic equations of this model are composed of two ordinary
@@ -16,23 +16,23 @@
         \dot{V} &= d \, \tau (-f V^3 + e V^2 + g V + \alpha W + \gamma I) \\
         \dot{W} &= \dfrac{d}{	au}\,\,(c V^2 + b V - \beta W + a)
         \end{align}
-```
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Other parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.
Citations: FitzHugh, R., Impulses and physiological states in theoretical models of nerve membrane, Biophysical Journal 1: 445, 1961.
Nagumo et.al, An Active Pulse Transmission Line Simulating Nerve Axon, Proceedings of the IRE 50: 2061, 1962.
Stefanescu, R., Jirsa, V.K. Reduced representations of heterogeneous mixed neural networks with synaptic coupling. Physical Review E, 83, 2011.
Jirsa VK, Stefanescu R.  Neural population modes capture biologically realistic large-scale network dynamics. Bulletin of Mathematical Biology, 2010.
Stefanescu, R., Jirsa, V.K. A low dimensional description of globally coupled heterogeneous neural networks of excitatory and inhibitory neurons. PLoS Computational Biology, 4(11), 2008).
source
Neuroblox.HarmonicOscillator — Type
HarmonicOscillator(name, namespace, ω, ζ, k, h)
+```
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Other parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.
Citations: FitzHugh, R., Impulses and physiological states in theoretical models of nerve membrane, Biophysical Journal 1: 445, 1961.
Nagumo et.al, An Active Pulse Transmission Line Simulating Nerve Axon, Proceedings of the IRE 50: 2061, 1962.
Stefanescu, R., Jirsa, V.K. Reduced representations of heterogeneous mixed neural networks with synaptic coupling. Physical Review E, 83, 2011.
Jirsa VK, Stefanescu R.  Neural population modes capture biologically realistic large-scale network dynamics. Bulletin of Mathematical Biology, 2010.
Stefanescu, R., Jirsa, V.K. A low dimensional description of globally coupled heterogeneous neural networks of excitatory and inhibitory neurons. PLoS Computational Biology, 4(11), 2008).
source
Neuroblox.HarmonicOscillator — Type
HarmonicOscillator(name, namespace, ω, ζ, k, h)
 
 Create a harmonic oscillator blox with the specified parameters.
 The formal definition of this blox is:
\[\frac{dx}{dt} = y-(2*\omega*\zeta*x)+ k*(2/\pi)*(atan((\sum{jcn})/h)
-\frac{dy}{dt} = -(\omega^2)*x\]
where ``jcn`` is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • ω: Base frequency. Note the default value is scaled to give oscillations in milliseconds to match other blocks.
  • ζ: Damping ratio.
  • k: Gain.
  • h: Threshold.
source
Neuroblox.JansenRit — Type
JansenRit(name, namespace, τ, H, λ, r, cortical, delayed)
+\frac{dy}{dt} = -(\omega^2)*x\]
where ``jcn`` is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • ω: Base frequency. Note the default value is scaled to give oscillations in milliseconds to match other blocks.
  • ζ: Damping ratio.
  • k: Gain.
  • h: Threshold.
source
Neuroblox.JansenRit — Type
JansenRit(name, namespace, τ, H, λ, r, cortical, delayed)
 
 Create a Jansen Rit blox as described in Liu et al.
 The formal definition of this blox is:
\[\frac{dx}{dt} = y-\frac{2}{\tau}x
-\frac{dy}{dt} = -\frac{x}{\tau^2} + \frac{H}{\tau} [\frac{2\lambda}{1+\text{exp}(-r*\sum{jcn})} - \lambda]\]
where $jcn$ is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • τ: Time constant. This is changed from the original source as the time constant was in seconds, while all our blocks are in milliseconds.
  • H: See equation for use.
  • λ: See equation for use.
  • r: See equation for use.
  • cortical: Boolean to determine whether to use cortical or subcortical parameters. Specifying any of the parameters above will override this.
  • delayed: Boolean to indicate whether states are delayed
Citations:
  1. Liu C, Zhou C, Wang J, Fietkiewicz C, Loparo KA. The role of coupling connections in a model of the cortico-basal ganglia-thalamocortical neural loop for the generation of beta oscillations. Neural Netw. 2020 Mar;123:381-392. doi: 10.1016/j.neunet.2019.12.021.
source
Neuroblox.LarterBreakspear — Type
LarterBreakspear(name, namespace, ...)
+\frac{dy}{dt} = -\frac{x}{\tau^2} + \frac{H}{\tau} [\frac{2\lambda}{1+\text{exp}(-r*\sum{jcn})} - \lambda]\]
where $jcn$ is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • τ: Time constant. This is changed from the original source as the time constant was in seconds, while all our blocks are in milliseconds.
  • H: See equation for use.
  • λ: See equation for use.
  • r: See equation for use.
  • cortical: Boolean to determine whether to use cortical or subcortical parameters. Specifying any of the parameters above will override this.
  • delayed: Boolean to indicate whether states are delayed
Citations:
  1. Liu C, Zhou C, Wang J, Fietkiewicz C, Loparo KA. The role of coupling connections in a model of the cortico-basal ganglia-thalamocortical neural loop for the generation of beta oscillations. Neural Netw. 2020 Mar;123:381-392. doi: 10.1016/j.neunet.2019.12.021.
source
Neuroblox.LarterBreakspear — Type
LarterBreakspear(name, namespace, ...)
 
 Create a Larter Breakspear blox described in Endo et al. For a full list of the parameters used see the reference.
-If you need to modify the parameters, see Chesebro et al. and van Nieuwenhuizen et al. for physiological ranges.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Other parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.
Citations:
  1. Endo H, Hiroe N, Yamashita O. Evaluation of Resting Spatio-Temporal Dynamics of a Neural Mass Model Using Resting fMRI Connectivity and EEG Microstates. Front Comput Neurosci. 2020 Jan 17;13:91. doi: 10.3389/fncom.2019.00091.
  2. Chesebro AG, Mujica-Parodi LR, Weistuch C. Ion gradient-driven bifurcations of a multi-scale neuronal model. Chaos Solitons Fractals. 2023 Feb;167:113120. doi: 10.1016/j.chaos.2023.113120. 
  3. van Nieuwenhuizen, H, Chesebro, AG, Polis, C, Clarke, K, Strey, HH, Weistuch, C, Mujica-Parodi, LR. Ketosis regulates K+ ion channels, strengthening brain-wide signaling disrupted by age. Preprint. bioRxiv 2023.05.10.540257; doi: https://doi.org/10.1101/2023.05.10.540257. 
source
Neuroblox.OUBlox — Type
Ornstein-Uhlenbeck process Blox
variables:     x(t):  value     jcn:   input  parameters:     τ:      relaxation time 	μ:      average value 	σ:      random noise (variance of OU process is τ*σ^2/2) returns:     an ODE System (but with brownian parameters)
source
Neuroblox.OUCouplingBlox — Type
Ornstein-Uhlenbeck Coupling Blox This blox takes an input and multiplies that input with a OU process of mean μ and variance τ*σ^2/2
This blox allows to create edges that have fluctuating weights
variables:     x(t):  value     jcn:   input  parameters:     τ:      relaxation time 	μ:      average value 	σ:      random noise (variance of OU process is τ*σ^2/2) returns:     an ODE System (but with brownian parameters)
source
Neuroblox.WilsonCowan — Type
WilsonCowan(name, namespace, τ_E, τ_I, a_E, a_I, c_EE, c_IE, c_EI, c_II, θ_E, θ_I, η)
+If you need to modify the parameters, see Chesebro et al. and van Nieuwenhuizen et al. for physiological ranges.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Other parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.
Citations:
  1. Endo H, Hiroe N, Yamashita O. Evaluation of Resting Spatio-Temporal Dynamics of a Neural Mass Model Using Resting fMRI Connectivity and EEG Microstates. Front Comput Neurosci. 2020 Jan 17;13:91. doi: 10.3389/fncom.2019.00091.
  2. Chesebro AG, Mujica-Parodi LR, Weistuch C. Ion gradient-driven bifurcations of a multi-scale neuronal model. Chaos Solitons Fractals. 2023 Feb;167:113120. doi: 10.1016/j.chaos.2023.113120. 
  3. van Nieuwenhuizen, H, Chesebro, AG, Polis, C, Clarke, K, Strey, HH, Weistuch, C, Mujica-Parodi, LR. Ketosis regulates K+ ion channels, strengthening brain-wide signaling disrupted by age. Preprint. bioRxiv 2023.05.10.540257; doi: https://doi.org/10.1101/2023.05.10.540257. 
source
Neuroblox.OUBlox — Type
Ornstein-Uhlenbeck process Blox
variables:     x(t):  value     jcn:   input  parameters:     τ:      relaxation time 	μ:      average value 	σ:      random noise (variance of OU process is τ*σ^2/2) returns:     an ODE System (but with brownian parameters)
source
Neuroblox.WilsonCowan — Type
WilsonCowan(name, namespace, τ_E, τ_I, a_E, a_I, c_EE, c_IE, c_EI, c_II, θ_E, θ_I, η)
 
 Create a standard Wilson Cowan blox.
 The formal definition of this blox is:
\[\frac{dE}{dt} = \frac{-E}{\tau_E} + \frac{1}{1 + \text{exp}(-a_E*(c_{EE}*E - c_{IE}*I - \theta_E + \eta*(\sum{jcn}))}
-\frac{dI}{dt} = \frac{-I}{\tau_I} + \frac{1}{1 + exp(-a_I*(c_{EI}*E - c_{II}*I - \theta_I)}\]
where $jcn$ is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Others: See equation for use.
source
Neuroblox.WinnerTakeAllBlox — Type
WinnerTakeAllBlox
Creates a winner-take-all local circuit found in neocortex, typically 5 pyramidal (excitatory) neurons send synapses to a single interneuron (inhibitory) and receive feedback inhibition from that interneuron.
source
LinearAlgebra.eigen — Method
function LinearAlgebra.eigen(M::Matrix{Dual{T, P, np}}) where {T, P, np}
+\frac{dI}{dt} = \frac{-I}{\tau_I} + \frac{1}{1 + exp(-a_I*(c_{EI}*E - c_{II}*I - \theta_I)}\]
where $jcn$ is any input to the blox.
Arguments:
  • name: Name given to ODESystem object within the blox.
  • namespace: Additional namespace above name if needed for inheritance.
  • Others: See equation for use.
source
Neuroblox.WinnerTakeAllBlox — Type
WinnerTakeAllBlox
Creates a winner-take-all local circuit found in neocortex, typically 5 pyramidal (excitatory) neurons send synapses to a single interneuron (inhibitory) and receive feedback inhibition from that interneuron.
source
LinearAlgebra.eigen — Method
function LinearAlgebra.eigen(M::Matrix{Dual{T, P, np}}) where {T, P, np}
 
 Dispatch of LinearAlgebra.eigen for dual matrices with complex numbers. Make the eigenvalue decomposition 
 amenable to automatic differentiation. To do so compute the analytical derivative of eigenvalues
@@ -42,7 +42,7 @@
 - `M`: matrix of type Dual of which to compute the eigenvalue decomposition. 
 
 Returns:
-- `Eigen(evals, evecs)`: eigenvalue decomposition returned as type LinearAlgebra.Eigen
source
Neuroblox.ARVTarget — Method
ARVTarget Time series data is bandpass filtered and then the power spectrum is computed for a given time interval (control bin), returned as the average value of the power spectral density within a certain frequency band ([lb, ub]).
source
Neuroblox.CDVTarget — Method
CDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Circular difference is quantified as the angle of circular_location.
source
Neuroblox.PDVTarget — Method
PDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Phase deviation is quantified as the angle difference between a given set of signals.
source
Neuroblox.PLVTarget — Method
PLVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians.
source
Neuroblox.addnontunableparams — Method
function addnontunableparams(param, model)
+- `Eigen(evals, evecs)`: eigenvalue decomposition returned as type LinearAlgebra.Eigen
source
Neuroblox.ARVTarget — Method
ARVTarget Time series data is bandpass filtered and then the power spectrum is computed for a given time interval (control bin), returned as the average value of the power spectral density within a certain frequency band ([lb, ub]).
source
Neuroblox.CDVTarget — Method
CDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Circular difference is quantified as the angle of circular_location.
source
Neuroblox.PDVTarget — Method
PDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Phase deviation is quantified as the angle difference between a given set of signals.
source
Neuroblox.PLVTarget — Method
PLVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians.
source
Neuroblox.addnontunableparams — Method
function addnontunableparams(param, model)
 
 Function adds parameters of a model that were not marked as tunable to a list of tunable parameters
 and respects the MTK ordering of parameters.
@@ -52,12 +52,12 @@
 - `sys`: MTK system
 
 Returns:
-- `completeparamlist`: complete parameter list of a system, including those that were not tagged as tunable
source
Neuroblox.bandpassfilter — Method
bandpassfilter takes in time series data and bandpass filters it. It has the following inputs:     data: time series data     lb: minimum cut-off frequency     ub: maximum cut-off frequency     fs: sampling frequency     order: filter order
source
Neuroblox.boldsignal — Method
Arguments:
  • name: Name given to ODESystem object within the blox.
  • lnϵ  : logarithm of ratio of intra- to extra-vascular signal
NB: the prefix ln of the variables ν, q as well as the parameters ϵ denotes their transformation into logarithmic space to enforce their positivity.
Citations:
  1. Stephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040
  2. Hofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407
source
Neuroblox.complexwavelet — Function
complexwavelet creates a complex morlet wavelet by windowing a complex sine wave with a Gaussian taper.  The morlet wavelet is a special case of a bandpass filter in which the frequency response is Gaussian-shaped. Convolution with a complex wavelet is equivalent to performing a Hilbert transform of a bandpass filtered signal.
It has the following inputs:     data: time series data      dt  : data sampling rate      lb  : lower bound wavelet frequency (in Hz)     ub  : upper bound wavelet frequency (in Hz)     a   : amplitude of the Gaussian taper, default is 1     n   : number of wavelet cycles of the Gaussian taper, defines the trade-off between temporal precision and frequency precision           larger n gives better frequency precision at the cost of temporal precision           default is 6 Hz     m   : x-axis offset, default is 0     num_wavelets : number of wavelets to create, default is 5
And outputs:     complex_wavelet : a family of complex morlet wavelets
source
Neuroblox.csd2mar — Method
This function converts a cross-spectral density (CSD) into a multivariate auto-regression (MAR) model. It first transforms the CSD into its cross-correlation function (Wiener-Kinchine theorem) and then computes the MAR model coefficients. csd       : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared) w         : frequencies dt        : time step size p         : number of time steps of auto-regressive model
This function returns coeff     : array of length p of coefficient matrices of size sqrt(N)xsqrt(N) noise_cov : noise covariance matrix
source
Neuroblox.csd_approx — Method
This function implements equation 2 of the spectral DCM paper, Friston et al. 2014 "A DCM for resting state fMRI".
+- `completeparamlist`: complete parameter list of a system, including those that were not tagged as tunable
source
Neuroblox.bandpassfilter — Method
bandpassfilter takes in time series data and bandpass filters it. It has the following inputs:     data: time series data     lb: minimum cut-off frequency     ub: maximum cut-off frequency     fs: sampling frequency     order: filter order
source
Neuroblox.boldsignal — Method
Arguments:
  • name: Name given to ODESystem object within the blox.
  • lnϵ  : logarithm of ratio of intra- to extra-vascular signal
NB: the prefix ln of the variables ν, q as well as the parameters ϵ denotes their transformation into logarithmic space to enforce their positivity.
Citations:
  1. Stephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040
  2. Hofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407
source
Neuroblox.complexwavelet — Function
complexwavelet creates a complex morlet wavelet by windowing a complex sine wave with a Gaussian taper.  The morlet wavelet is a special case of a bandpass filter in which the frequency response is Gaussian-shaped. Convolution with a complex wavelet is equivalent to performing a Hilbert transform of a bandpass filtered signal.
It has the following inputs:     data: time series data      dt  : data sampling rate      lb  : lower bound wavelet frequency (in Hz)     ub  : upper bound wavelet frequency (in Hz)     a   : amplitude of the Gaussian taper, default is 1     n   : number of wavelet cycles of the Gaussian taper, defines the trade-off between temporal precision and frequency precision           larger n gives better frequency precision at the cost of temporal precision           default is 6 Hz     m   : x-axis offset, default is 0     num_wavelets : number of wavelets to create, default is 5
And outputs:     complex_wavelet : a family of complex morlet wavelets
source
Neuroblox.csd2mar — Method
This function converts a cross-spectral density (CSD) into a multivariate auto-regression (MAR) model. It first transforms the CSD into its cross-correlation function (Wiener-Kinchine theorem) and then computes the MAR model coefficients. csd       : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared) w         : frequencies dt        : time step size p         : number of time steps of auto-regressive model
This function returns coeff     : array of length p of coefficient matrices of size sqrt(N)xsqrt(N) noise_cov : noise covariance matrix
source
Neuroblox.csd_approx — Method
This function implements equation 2 of the spectral DCM paper, Friston et al. 2014 "A DCM for resting state fMRI".
 Note that nomenclature is taken from SPM12 code and it does not seem to coincide with the spectral DCM paper's nomenclature. 
 For instance, Gu should represent the spectral component due to external input according to the paper. However, in the code this represents
 the hidden state fluctuations (which are called Gν in the paper).
 Gn in the code corresponds to Ge in the paper, i.e. the observation noise. In the code global and local components are defined, no such distinction
-is discussed in the paper. In fact the parameter γ, corresponding to local component is not present in the paper.
source
Neuroblox.get_dynamic_states — Method
function get_dynamic_states(sys)
+is discussed in the paper. In fact the parameter γ, corresponding to local component is not present in the paper.
source
Neuroblox.get_dynamic_states — Method
function get_dynamic_states(sys)
 
 Function extracts states from the system that are dynamic variables, 
 get also indices of external inputs (u(t)) and measurements (like bold(t))
@@ -65,11 +65,10 @@
 - `sys`: MTK system
 
 Returns:
-- `sts`  : states of the system that are neither external inputs nor measurements, i.e. these are the dynamic states
-- `idx_u`: indices of states that represent external inputs
-- `idx_m`: indices of states that represent measurements
source
Neuroblox.idft — Method
Plain implementation of idft because AD dispatch versions for ifft don't work still!
source
Neuroblox.inner_namespaceof — Method
Returns the complete namespace EXCLUDING the outermost (highest) level.
+- `sts`: states/unknowns of the system that are neither external inputs nor measurements, i.e. these are the dynamic states
+- `idx`: indices of these states
source
Neuroblox.idft — Method
Plain implementation of idft because AD dispatch versions for ifft don't work still!
source
Neuroblox.inner_namespaceof — Method
Returns the complete namespace EXCLUDING the outermost (highest) level.
 This is useful for manually preparing equations (e.g. connections, see BloxConnector),
-that will later be composed and will automatically get the outermost namespace.
source
Neuroblox.input_equations — Method
Returns the equations for all input variables of a system, 
+that will later be composed and will automatically get the outermost namespace.
source
Neuroblox.input_equations — Method
Returns the equations for all input variables of a system, 
 assuming they have a form like : `sys.input_variable ~ ...`
 so only the input appears on the LHS.
 
@@ -79,18 +78,18 @@
 
 If blox isa AbstractComponent, it is assumed that it contains a `connector` field,
 which holds a `BloxConnector` object with all relevant connections 
-from lower levels and this level.
source
Neuroblox.learningrate — Method
This function computes learning rate. It has the following inputs:     outcomes: vector of 1's and 0's for behavioral outcomes     windows: number of windows to split the outcome data into And the following outputs:     rate: the learning rate across each window
source
Neuroblox.mar2csd — Method
This function converts multivariate auto-regression (MAR) model parameters to a cross-spectral density (CSD). A     : coefficients of MAR model, array of length p, each element contains the regression coefficients for that particular time-lag. Σ     : noise covariance matrix of MAR p     : number of time lags freqs : frequencies at which to evaluate the CSD sf    : sampling frequency
This function returns: csd   : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared)
source
Neuroblox.mar_ml — Method
Maximum likelihood estimator of a multivariate, or vector auto-regressive model.     y : MxN Data matrix where M is number of samples and N is number of dimensions     p : time lag parameter, also called order of MAR model     return values     mar["A"] : model parameters is a NxNxP tensor, i.e. one NxN parameter matrix for each time bin k ∈ {1,...,p}     mar["Σ"] : noise covariance matrix
source
Neuroblox.learningrate — Method
This function computes learning rate. It has the following inputs:     outcomes: vector of 1's and 0's for behavioral outcomes     windows: number of windows to split the outcome data into And the following outputs:     rate: the learning rate across each window
source
Neuroblox.mar2csd — Method
This function converts multivariate auto-regression (MAR) model parameters to a cross-spectral density (CSD). A     : coefficients of MAR model, array of length p, each element contains the regression coefficients for that particular time-lag. Σ     : noise covariance matrix of MAR p     : number of time lags freqs : frequencies at which to evaluate the CSD sf    : sampling frequency
This function returns: csd   : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared)
source
Neuroblox.mar_ml — Method
Maximum likelihood estimator of a multivariate, or vector auto-regressive model.     y : MxN Data matrix where M is number of samples and N is number of dimensions     p : time lag parameter, also called order of MAR model     return values     mar["A"] : model parameters is a NxNxP tensor, i.e. one NxN parameter matrix for each time bin k ∈ {1,...,p}     mar["Σ"] : noise covariance matrix
source
Neuroblox.matlab_norm — Method
function matlab_norm(A, p)
 
 Simple helper function to implement the norm of a matrix that is equivalent to the one given in MATLAB for order=1, 2, Inf. 
 This is needed for the reproduction of the exact same results of SPM12.
 
 Arguments:
 - `A`: matrix
-- `p`: order of norm
source
Neuroblox.paramscoping — Method
function paramscoping(;kwargs...)
 
 Scope arguments that are already a symbolic model parameter thereby keep the correct namespace 
 and make those that are not yet symbolic a symbol.
-Keyword arguments are used, because parameter definition require names, not just values.
source
Neuroblox.phase_cos_blox — Method
phasecosblox is creating a cos with angular frequency ω and variable phase phaseinter has the following parameters:     ω: angular frequency     t: time     phaseinter: a function that returns phase as a function of time and returns:     the resulting value
Usage:     phaseint = phaseinter(0:0.1:50,phasedata)     phaseout(t) = phasecosblox(0.1,t,phaseint)     which is now a function of time and can be used in an input blox     you can also use the dot operator to calculate time-series     signal = phaseout.(collect(0:0.01:50))
source
Neuroblox.phase_inter — Method
phaseinter is creating a function that interpolates the phase data for any time given phaseinter has the following parameters:     phaserange:  a range, e.g. 0:0.1:50 which should reflect the time points of the data     phasedata: phase at equidistant time points and returns:     an function that returns an interpolated phase for t in range
source
Neuroblox.phase_sin_blox — Method
phasesinblox is creating a sin with angular frequency ω and variable phase phaseinter has the following parameters:     ω: angular frequency     t: time     phaseinter: a function that returns phase as a function of time and returns:     the resulting value
Usage:     phaseint = phaseinter(0:0.1:50,phasedata)     phaseout(t) = phasesinblox(0.1,t,phaseint)     which is now a function of time and can be used in an input blox     you can also use the dot operator to calculate time-series     signal = phaseout.(collect(0:0.01:50))
source
Neuroblox.random_initials — Method
random_initials creates a vector of random initial conditions for an ODESystem that is composed of a list of blox.  The function finds the initial conditions in the blox and then sets a random value in between range tuple given for that state.
It has the following inputs:     odesys: ODESystem     blox  : list of blox
And outputs:     u0 : Float64 vector of initial conditions
source
Neuroblox.setup_sDCM — Method
function setup_sDCM(data, stateevolutionmodel, initcond, csdsetup, priors, hyperpriors, indices)
+Keyword arguments are used, because parameter definition require names, not just values.
source
Neuroblox.phase_cos_blox — Method
phasecosblox is creating a cos with angular frequency ω and variable phase phaseinter has the following parameters:     ω: angular frequency     t: time     phaseinter: a function that returns phase as a function of time and returns:     the resulting value
Usage:     phaseint = phaseinter(0:0.1:50,phasedata)     phaseout(t) = phasecosblox(0.1,t,phaseint)     which is now a function of time and can be used in an input blox     you can also use the dot operator to calculate time-series     signal = phaseout.(collect(0:0.01:50))
source
Neuroblox.phase_inter — Method
phaseinter is creating a function that interpolates the phase data for any time given phaseinter has the following parameters:     phaserange:  a range, e.g. 0:0.1:50 which should reflect the time points of the data     phasedata: phase at equidistant time points and returns:     an function that returns an interpolated phase for t in range
source
Neuroblox.phase_sin_blox — Method
phasesinblox is creating a sin with angular frequency ω and variable phase phaseinter has the following parameters:     ω: angular frequency     t: time     phaseinter: a function that returns phase as a function of time and returns:     the resulting value
Usage:     phaseint = phaseinter(0:0.1:50,phasedata)     phaseout(t) = phasesinblox(0.1,t,phaseint)     which is now a function of time and can be used in an input blox     you can also use the dot operator to calculate time-series     signal = phaseout.(collect(0:0.01:50))
source
Neuroblox.random_initials — Method
random_initials creates a vector of random initial conditions for an ODESystem that is composed of a list of blox.  The function finds the initial conditions in the blox and then sets a random value in between range tuple given for that state.
It has the following inputs:     odesys: ODESystem     blox  : list of blox
And outputs:     u0 : Float64 vector of initial conditions
source
Neuroblox.setup_sDCM — Method
function setup_sDCM(data, stateevolutionmodel, initcond, csdsetup, priors, hyperpriors, indices)
 
 Interface function to performs variational inference to fit model parameters to empirical cross spectral density.
 The current implementation provides a Variational Laplace fit (see function above `variationalbayes`).
@@ -110,14 +109,14 @@
 - `hyperpriors` : dataframe of parameters with the following columns:
 -- `Πλ_pr`      : prior precision matrix for λ hyperparameter(s)
 -- `μλ_pr`      : prior mean(s) for λ hyperparameter(s)
-- `indices`  : indices to separate model parameters from other parameters. Needed for the computation of AD gradient.
source
Neuroblox.spm_logdet — Method
function spm_logdet(M)
+- `indices`  : indices to separate model parameters from other parameters. Needed for the computation of AD gradient.
source
Neuroblox.spm_logdet — Method
function spm_logdet(M)
 
 SPM12 style implementation of the logarithm of the determinant of a matrix.
 
 Arguments:
-- `M`: matrix
source
Neuroblox.vecparam — Method
vecparam(param::OrderedDict)
 
 Function to flatten an ordered dictionary of model parameters and return a simple list of parameter values.
 
 Arguments:
-- `param`: dictionary of model parameters (may contain numbers and lists of numbers)
source

+- `param`: dictionary of model parameters (may contain numbers and lists of numbers)
source
diff --git a/dev/assets/Manifest.toml b/dev/assets/Manifest.toml index b63e983..a7713f2 100644 --- a/dev/assets/Manifest.toml +++ b/dev/assets/Manifest.toml @@ -305,9 +305,9 @@ version = "5.9.0" ODEInterface = "54ca160b-1b9f-5127-a996-1867f4bc2a2c" [[deps.BufferedStreams]] -git-tree-sha1 = "4ae47f9a4b1dc19897d3743ff13685925c5202ec" +git-tree-sha1 = "6863c5b7fc997eadcabdbaf6c5f201dc30032643" uuid = "e1450e63-4bb3-523b-b2a4-4ffa8c0fd77d" -version = "1.2.1" +version = "1.2.2" [[deps.Bzip2_jll]] deps = ["Artifacts", "JLLWrappers", "Libdl", "Pkg"] @@ -2036,9 +2036,9 @@ version = "9.26.0" [[deps.ModelingToolkitStandardLibrary]] deps = ["ChainRulesCore", "DiffEqBase", "IfElse", "LinearAlgebra", "ModelingToolkit", "Symbolics"] -git-tree-sha1 = "862e6f77409c3ee33ade20b0242b10225fcc5e5a" +git-tree-sha1 = "53c2eb67a861a4c2bb5d37b3526659e033dc5c31" uuid = "16a59e39-deab-5bd0-87e4-056b12336739" -version = "2.9.0" +version = "2.10.0" [[deps.MozillaCACerts_jll]] uuid = "14a3606d-f60d-562e-9121-12d972cd8159" @@ -2313,9 +2313,9 @@ version = "1.6.3" [[deps.OrdinaryDiffEq]] deps = ["ADTypes", "Adapt", "ArrayInterface", "DataStructures", "DiffEqBase", "DocStringExtensions", "EnumX", "ExponentialUtilities", "FastBroadcast", "FastClosures", "FillArrays", "FiniteDiff", "ForwardDiff", "FunctionWrappersWrappers", "IfElse", "InteractiveUtils", "LineSearches", "LinearAlgebra", "LinearSolve", "Logging", "MacroTools", "MuladdMacro", "NonlinearSolve", "Polyester", "PreallocationTools", "PrecompileTools", "Preferences", "RecursiveArrayTools", "Reexport", "SciMLBase", "SciMLOperators", "SciMLStructures", "SimpleNonlinearSolve", "SimpleUnPack", "SparseArrays", "SparseDiffTools", "Static", "StaticArrayInterface", "StaticArrays", "TruncatedStacktraces"] -git-tree-sha1 = "0c5d89483f9538efedb3f1c1b72e14d5f65830b0" +git-tree-sha1 = "a8b2d333cd90562b58b977b4033739360b37fb1f" uuid = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" -version = "6.86.0" +version = "6.87.0" [[deps.PCRE2_jll]] deps = ["Artifacts", "Libdl"] @@ -2721,9 +2721,9 @@ version = "0.1.0" [[deps.SciMLBase]] deps = ["ADTypes", "Accessors", "ArrayInterface", "CommonSolve", "ConstructionBase", "Distributed", "DocStringExtensions", "EnumX", "Expronicon", "FunctionWrappersWrappers", "IteratorInterfaceExtensions", "LinearAlgebra", "Logging", "Markdown", "PrecompileTools", "Preferences", "Printf", "RecipesBase", "RecursiveArrayTools", "Reexport", "RuntimeGeneratedFunctions", "SciMLOperators", "SciMLStructures", "StaticArraysCore", "Statistics", "SymbolicIndexingInterface", "Tables"] -git-tree-sha1 = "8501bb749cc6684c86b6f113e17f2e6d60ccf11c" +git-tree-sha1 = "380a059a9fd18a56d98e50ed98d40e1c1202e662" uuid = "0bca4576-84f4-4d90-8ffe-ffa030f20462" -version = "2.45.1" +version = "2.46.0" [deps.SciMLBase.extensions] SciMLBaseChainRulesCoreExt = "ChainRulesCore" @@ -3024,9 +3024,9 @@ version = "2.2.0" [[deps.StochasticDiffEq]] deps = ["Adapt", "ArrayInterface", "DataStructures", "DiffEqBase", "DiffEqNoiseProcess", "DocStringExtensions", "FiniteDiff", "ForwardDiff", "JumpProcesses", "LevyArea", "LinearAlgebra", "Logging", "MuladdMacro", "NLsolve", "OrdinaryDiffEq", "Random", "RandomNumbers", "RecursiveArrayTools", "Reexport", "SciMLBase", "SciMLOperators", "SparseArrays", "SparseDiffTools", "StaticArrays", "UnPack"] -git-tree-sha1 = "5237f2ebdf3b6b7ad2ec635440b59a390988feeb" +git-tree-sha1 = "b47f8ccc5bd06d5f7a643bf6671365ab9d6595d9" uuid = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0" -version = "6.66.0" +version = "6.67.0" [[deps.StrideArraysCore]] deps = ["ArrayInterface", "CloseOpenIntervals", "IfElse", "LayoutPointers", "LinearAlgebra", "ManualMemory", "SIMDTypes", "Static", "StaticArrayInterface", "ThreadingUtilities"] @@ -3232,9 +3232,9 @@ uuid = "d5829a12-d9aa-46ab-831f-fb7c9ab06edf" version = "0.2.1" [[deps.Tricks]] -git-tree-sha1 = "eae1bb484cd63b36999ee58be2de6c178105112f" +git-tree-sha1 = "7822b97e99a1672bfb1b49b668a6d46d58d8cbcb" uuid = "410a4b4d-49e4-4fbc-ab6d-cb71b17b3775" -version = "0.1.8" +version = "0.1.9" [[deps.TruncatedStacktraces]] deps = ["InteractiveUtils", "MacroTools", "Preferences"] diff --git a/dev/getting_started/8e049f77.svg b/dev/getting_started/682067d5.svg similarity index 96% rename from dev/getting_started/8e049f77.svg rename to dev/getting_started/682067d5.svg index e936ec5..1ecbe23 100644 --- a/dev/getting_started/8e049f77.svg +++ b/dev/getting_started/682067d5.svg @@ -1,76 +1,76 @@ - + - + - + - + - + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/dev/getting_started/28da62dd.svg b/dev/getting_started/91fd99fd.svg similarity index 94% rename from dev/getting_started/28da62dd.svg rename to dev/getting_started/91fd99fd.svg index eb91595..7affd59 100644 --- a/dev/getting_started/28da62dd.svg +++ b/dev/getting_started/91fd99fd.svg @@ -1,54 +1,54 @@ - + - + - + - + - + - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/dev/getting_started/index.html b/dev/getting_started/index.html index f065b8e..9892f59 100644 --- a/dev/getting_started/index.html +++ b/dev/getting_started/index.html @@ -22,7 +22,7 @@ \end{align} \]

To solve the system, we first create an ODEProblem and then solve it over the tspan of (0,100) using a stiff solver. The solution is saved every 0.1ms. The unit of time in Neuroblox is 1ms.

prob = ODEProblem(sys, [], (0.0, 100), [])
 sol = solve(prob, Rodas4(), saveat=0.1)
-plot(sol)
Example block output

Example 2 : Building a Brain Circuit from literature using Neural Mass Models

In this example, we will construct a Parkinsons model from eight Jansen-Rit Neural Mass Models as described in Liu et al. (2020). DOI: 10.1016/j.neunet.2019.12.021. The Jansen-Rit Neural Mass model is defined by the following differential equations:

\[\frac{dx}{dt} = y-\frac{2}{\tau}x +plot(sol)Example block output

Example 2 : Building a Brain Circuit from literature using Neural Mass Models

In this example, we will construct a Parkinsons model from eight Jansen-Rit Neural Mass Models as described in Liu et al. (2020). DOI: 10.1016/j.neunet.2019.12.021. The Jansen-Rit Neural Mass model is defined by the following differential equations:

\[\frac{dx}{dt} = y-\frac{2}{\tau}x \frac{dy}{dt} = -\frac{x}{\tau^2} + \frac{H}{\tau} [\frac{2\lambda}{1+\text{exp}(-r*\sum{jcn})} - \lambda]\]

using Neuroblox
 using DifferentialEquations
 using Graphs
@@ -121,4 +121,4 @@
  1.0
  1.0

We select an algorithm and solve the system

alg = Tsit5()
 sol_dde_no_delays = solve(prob, alg, saveat=1)
-plot(sol_dde_no_delays)
Example block output

In a later tutorial, we will show how to introduce edge delays.

+plot(sol_dde_no_delays)Example block output

In a later tutorial, we will show how to introduce edge delays.

diff --git a/dev/index.html b/dev/index.html index 0ba2667..8acfb63 100644 --- a/dev/index.html +++ b/dev/index.html @@ -4,4 +4,4 @@ using PkgAuthentication PkgAuthentication.install("juliahub.com") Pkg.Registry.add() -Pkg.add("Neuroblox")

Licensing

Neuroblox is free for non-commerical and academic use. For full details of the license, please see the Neuroblox EULA. For commercial use, get in contact with sales@neuroblox.org.

+Pkg.add("Neuroblox")

Licensing

Neuroblox is free for non-commerical and academic use. For full details of the license, please see the Neuroblox EULA. For commercial use, get in contact with sales@neuroblox.org.

diff --git a/dev/objects.inv b/dev/objects.inv index 2d27040..3011d5e 100644 Binary files a/dev/objects.inv and b/dev/objects.inv differ diff --git a/dev/release_notes/index.html b/dev/release_notes/index.html index 8ef402b..496cae5 100644 --- a/dev/release_notes/index.html +++ b/dev/release_notes/index.html @@ -1,2 +1,2 @@ -Release Notes · Neuroblox
GitHub
+Release Notes · Neuroblox
GitHub
diff --git a/dev/search_index.js b/dev/search_index.js index b713121..04acb8e 100644 --- a/dev/search_index.js +++ b/dev/search_index.js @@ -1,3 +1,3 @@ var documenterSearchIndex = {"docs": -[{"location":"release_notes/#Release-Notes","page":"Release Notes","title":"Release Notes","text":"","category":"section"},{"location":"release_notes/#v0.3","page":"Release Notes","title":"v0.3","text":"","category":"section"},{"location":"release_notes/","page":"Release Notes","title":"Release Notes","text":"Initial Release!","category":"page"},{"location":"api/#API-Documentation","page":"API","title":"API Documentation","text":"","category":"section"},{"location":"api/","page":"API","title":"API","text":"Modules = [Neuroblox]","category":"page"},{"location":"api/#Neuroblox.BalloonModel","page":"API","title":"Neuroblox.BalloonModel","text":"Arguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nlnκ: logarithmic prefactor to signal decay H[1], set to 0 for standard parameter value.\nlnτ: logarithmic prefactor to transit time H[3], set to 0 for standard parameter value.\nlnϵ: logarithm of ratio of intra- to extra-vascular signal\n\nNB: the prefix ln of the variables u, ν, q as well as the parameters κ, τ denotes their transformation into logarithmic space to enforce their positivity. This transformation is considered in the derivates of the model equations below. \n\nCitations:\n\nStephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040\nHofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.Generic2dOscillator","page":"API","title":"Neuroblox.Generic2dOscillator","text":"Generic2dOscillator(name, namespace, ...)\n\nThe Generic2dOscillator model is a generic dynamic system with two state\nvariables. The dynamic equations of this model are composed of two ordinary\ndifferential equations comprising two nullclines. The first nullcline is a\ncubic function as it is found in most neuron and population models; the\nsecond nullcline is arbitrarily configurable as a polynomial function up to\nsecond order. The manipulation of the latter nullcline's parameters allows\nto generate a wide range of different behaviours.\n\nEquations:\n\n```math\n \\begin{align}\n \\dot{V} &= d \\, \\tau (-f V^3 + e V^2 + g V + \\alpha W + \\gamma I) \\\\\n \\dot{W} &= \\dfrac{d}{\tau}\\,\\,(c V^2 + b V - \\beta W + a)\n \\end{align}\n```\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOther parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.\n\nCitations: FitzHugh, R., Impulses and physiological states in theoretical models of nerve membrane, Biophysical Journal 1: 445, 1961.\n\nNagumo et.al, An Active Pulse Transmission Line Simulating Nerve Axon, Proceedings of the IRE 50: 2061, 1962.\n\nStefanescu, R., Jirsa, V.K. Reduced representations of heterogeneous mixed neural networks with synaptic coupling. Physical Review E, 83, 2011.\n\nJirsa VK, Stefanescu R. Neural population modes capture biologically realistic large-scale network dynamics. Bulletin of Mathematical Biology, 2010.\n\nStefanescu, R., Jirsa, V.K. A low dimensional description of globally coupled heterogeneous neural networks of excitatory and inhibitory neurons. PLoS Computational Biology, 4(11), 2008).\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.HarmonicOscillator","page":"API","title":"Neuroblox.HarmonicOscillator","text":"HarmonicOscillator(name, namespace, ω, ζ, k, h)\n\nCreate a harmonic oscillator blox with the specified parameters.\nThe formal definition of this blox is:\n\nfracdxdt = y-(2*omega*zeta*x)+ k*(2pi)*(atan((sumjcn)h)\nfracdydt = -(omega^2)*x\n\nwhere ``jcn`` is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nω: Base frequency. Note the default value is scaled to give oscillations in milliseconds to match other blocks.\nζ: Damping ratio.\nk: Gain.\nh: Threshold.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.JansenRit","page":"API","title":"Neuroblox.JansenRit","text":"JansenRit(name, namespace, τ, H, λ, r, cortical, delayed)\n\nCreate a Jansen Rit blox as described in Liu et al.\nThe formal definition of this blox is:\n\nfracdxdt = y-frac2taux\nfracdydt = -fracxtau^2 + fracHtau frac2lambda1+textexp(-r*sumjcn) - lambda\n\nwhere jcn is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nτ: Time constant. This is changed from the original source as the time constant was in seconds, while all our blocks are in milliseconds.\nH: See equation for use.\nλ: See equation for use.\nr: See equation for use.\ncortical: Boolean to determine whether to use cortical or subcortical parameters. Specifying any of the parameters above will override this.\ndelayed: Boolean to indicate whether states are delayed\n\nCitations:\n\nLiu C, Zhou C, Wang J, Fietkiewicz C, Loparo KA. The role of coupling connections in a model of the cortico-basal ganglia-thalamocortical neural loop for the generation of beta oscillations. Neural Netw. 2020 Mar;123:381-392. doi: 10.1016/j.neunet.2019.12.021.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.JansenRitSPM12","page":"API","title":"Neuroblox.JansenRitSPM12","text":"Jansen-Rit model block for canonical micro circuit, analogous to the implementation in SPM12\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.LarterBreakspear","page":"API","title":"Neuroblox.LarterBreakspear","text":"LarterBreakspear(name, namespace, ...)\n\nCreate a Larter Breakspear blox described in Endo et al. For a full list of the parameters used see the reference.\nIf you need to modify the parameters, see Chesebro et al. and van Nieuwenhuizen et al. for physiological ranges.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOther parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.\n\nCitations:\n\nEndo H, Hiroe N, Yamashita O. Evaluation of Resting Spatio-Temporal Dynamics of a Neural Mass Model Using Resting fMRI Connectivity and EEG Microstates. Front Comput Neurosci. 2020 Jan 17;13:91. doi: 10.3389/fncom.2019.00091.\nChesebro AG, Mujica-Parodi LR, Weistuch C. Ion gradient-driven bifurcations of a multi-scale neuronal model. Chaos Solitons Fractals. 2023 Feb;167:113120. doi: 10.1016/j.chaos.2023.113120. \nvan Nieuwenhuizen, H, Chesebro, AG, Polis, C, Clarke, K, Strey, HH, Weistuch, C, Mujica-Parodi, LR. Ketosis regulates K+ ion channels, strengthening brain-wide signaling disrupted by age. Preprint. bioRxiv 2023.05.10.540257; doi: https://doi.org/10.1101/2023.05.10.540257. \n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.OUBlox","page":"API","title":"Neuroblox.OUBlox","text":"Ornstein-Uhlenbeck process Blox\n\nvariables: x(t): value jcn: input parameters: τ: relaxation time \tμ: average value \tσ: random noise (variance of OU process is τ*σ^2/2) returns: an ODE System (but with brownian parameters)\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.OUCouplingBlox","page":"API","title":"Neuroblox.OUCouplingBlox","text":"Ornstein-Uhlenbeck Coupling Blox This blox takes an input and multiplies that input with a OU process of mean μ and variance τ*σ^2/2\n\nThis blox allows to create edges that have fluctuating weights\n\nvariables: x(t): value jcn: input parameters: τ: relaxation time \tμ: average value \tσ: random noise (variance of OU process is τ*σ^2/2) returns: an ODE System (but with brownian parameters)\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.Striatum","page":"API","title":"Neuroblox.Striatum","text":"Subcortical blox\nall subcprtical blox used in cortico-striatal model are defined here\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.WilsonCowan","page":"API","title":"Neuroblox.WilsonCowan","text":"WilsonCowan(name, namespace, τ_E, τ_I, a_E, a_I, c_EE, c_IE, c_EI, c_II, θ_E, θ_I, η)\n\nCreate a standard Wilson Cowan blox.\nThe formal definition of this blox is:\n\nfracdEdt = frac-Etau_E + frac11 + textexp(-a_E*(c_EE*E - c_IE*I - theta_E + eta*(sumjcn))\nfracdIdt = frac-Itau_I + frac11 + exp(-a_I*(c_EI*E - c_II*I - theta_I)\n\nwhere jcn is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOthers: See equation for use.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.WinnerTakeAllBlox","page":"API","title":"Neuroblox.WinnerTakeAllBlox","text":"WinnerTakeAllBlox\n\nCreates a winner-take-all local circuit found in neocortex, typically 5 pyramidal (excitatory) neurons send synapses to a single interneuron (inhibitory) and receive feedback inhibition from that interneuron.\n\n\n\n\n\n","category":"type"},{"location":"api/#LinearAlgebra.eigen-Union{Tuple{Array{ForwardDiff.Dual{T, P, np}, 2}}, Tuple{np}, Tuple{P}, Tuple{T}} where {T, P, np}","page":"API","title":"LinearAlgebra.eigen","text":"function LinearAlgebra.eigen(M::Matrix{Dual{T, P, np}}) where {T, P, np}\n\nDispatch of LinearAlgebra.eigen for dual matrices with complex numbers. Make the eigenvalue decomposition \namenable to automatic differentiation. To do so compute the analytical derivative of eigenvalues\nand eigenvectors. \n\nArguments:\n- `M`: matrix of type Dual of which to compute the eigenvalue decomposition. \n\nReturns:\n- `Eigen(evals, evecs)`: eigenvalue decomposition returned as type LinearAlgebra.Eigen\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.ARVTarget-NTuple{6, Any}","page":"API","title":"Neuroblox.ARVTarget","text":"ARVTarget Time series data is bandpass filtered and then the power spectrum is computed for a given time interval (control bin), returned as the average value of the power spectral density within a certain frequency band ([lb, ub]).\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.CDVTarget-NTuple{5, Any}","page":"API","title":"Neuroblox.CDVTarget","text":"CDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Circular difference is quantified as the angle of circular_location.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.ControlError-NTuple{8, Any}","page":"API","title":"Neuroblox.ControlError","text":"ControlError Returns the control error (deviation of the actual value from the target value).\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.PDVTarget-NTuple{5, Any}","page":"API","title":"Neuroblox.PDVTarget","text":"PDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Phase deviation is quantified as the angle difference between a given set of signals.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.PLVTarget-NTuple{6, Any}","page":"API","title":"Neuroblox.PLVTarget","text":"PLVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.addnontunableparams-Tuple{Any, Any}","page":"API","title":"Neuroblox.addnontunableparams","text":"function addnontunableparams(param, model)\n\nFunction adds parameters of a model that were not marked as tunable to a list of tunable parameters\nand respects the MTK ordering of parameters.\n\nArguments:\n- `paramlist`: parameters of an MTK system that were tagged as tunable\n- `sys`: MTK system\n\nReturns:\n- `completeparamlist`: complete parameter list of a system, including those that were not tagged as tunable\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.bandpassfilter-Tuple{}","page":"API","title":"Neuroblox.bandpassfilter","text":"bandpassfilter takes in time series data and bandpass filters it. It has the following inputs: data: time series data lb: minimum cut-off frequency ub: maximum cut-off frequency fs: sampling frequency order: filter order\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.boldsignal-Tuple{}","page":"API","title":"Neuroblox.boldsignal","text":"Arguments:\n\nname: Name given to ODESystem object within the blox.\nlnϵ : logarithm of ratio of intra- to extra-vascular signal\n\nNB: the prefix ln of the variables ν, q as well as the parameters ϵ denotes their transformation into logarithmic space to enforce their positivity.\n\nCitations:\n\nStephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040\nHofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.complexwavelet","page":"API","title":"Neuroblox.complexwavelet","text":"complexwavelet creates a complex morlet wavelet by windowing a complex sine wave with a Gaussian taper. The morlet wavelet is a special case of a bandpass filter in which the frequency response is Gaussian-shaped. Convolution with a complex wavelet is equivalent to performing a Hilbert transform of a bandpass filtered signal.\n\nIt has the following inputs: data: time series data dt : data sampling rate lb : lower bound wavelet frequency (in Hz) ub : upper bound wavelet frequency (in Hz) a : amplitude of the Gaussian taper, default is 1 n : number of wavelet cycles of the Gaussian taper, defines the trade-off between temporal precision and frequency precision larger n gives better frequency precision at the cost of temporal precision default is 6 Hz m : x-axis offset, default is 0 num_wavelets : number of wavelets to create, default is 5\n\nAnd outputs: complex_wavelet : a family of complex morlet wavelets\n\n\n\n\n\n","category":"function"},{"location":"api/#Neuroblox.csd2mar-NTuple{4, Any}","page":"API","title":"Neuroblox.csd2mar","text":"This function converts a cross-spectral density (CSD) into a multivariate auto-regression (MAR) model. It first transforms the CSD into its cross-correlation function (Wiener-Kinchine theorem) and then computes the MAR model coefficients. csd : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared) w : frequencies dt : time step size p : number of time steps of auto-regressive model\n\nThis function returns coeff : array of length p of coefficient matrices of size sqrt(N)xsqrt(N) noise_cov : noise covariance matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.csd_approx-NTuple{4, Any}","page":"API","title":"Neuroblox.csd_approx","text":"This function implements equation 2 of the spectral DCM paper, Friston et al. 2014 \"A DCM for resting state fMRI\".\nNote that nomenclature is taken from SPM12 code and it does not seem to coincide with the spectral DCM paper's nomenclature. \nFor instance, Gu should represent the spectral component due to external input according to the paper. However, in the code this represents\nthe hidden state fluctuations (which are called Gν in the paper).\nGn in the code corresponds to Ge in the paper, i.e. the observation noise. In the code global and local components are defined, no such distinction\nis discussed in the paper. In fact the parameter γ, corresponding to local component is not present in the paper.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.get_dynamic_states-Tuple{Any}","page":"API","title":"Neuroblox.get_dynamic_states","text":"function get_dynamic_states(sys)\n\nFunction extracts states from the system that are dynamic variables, \nget also indices of external inputs (u(t)) and measurements (like bold(t))\nArguments:\n- `sys`: MTK system\n\nReturns:\n- `sts` : states of the system that are neither external inputs nor measurements, i.e. these are the dynamic states\n- `idx_u`: indices of states that represent external inputs\n- `idx_m`: indices of states that represent measurements\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.idft-Tuple{AbstractArray}","page":"API","title":"Neuroblox.idft","text":"Plain implementation of idft because AD dispatch versions for ifft don't work still!\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.inner_namespaceof-Tuple{Any}","page":"API","title":"Neuroblox.inner_namespaceof","text":"Returns the complete namespace EXCLUDING the outermost (highest) level.\nThis is useful for manually preparing equations (e.g. connections, see BloxConnector),\nthat will later be composed and will automatically get the outermost namespace.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.input_equations-Tuple{Any}","page":"API","title":"Neuroblox.input_equations","text":"Returns the equations for all input variables of a system, \nassuming they have a form like : `sys.input_variable ~ ...`\nso only the input appears on the LHS.\n\nInput equations are namespaced by the inner namespace of blox\nand then they are returned. This way during system `compose` downstream,\nthe higher-level namespaces will be added to them.\n\nIf blox isa AbstractComponent, it is assumed that it contains a `connector` field,\nwhich holds a `BloxConnector` object with all relevant connections \nfrom lower levels and this level.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.learningrate-Tuple{Any, Any}","page":"API","title":"Neuroblox.learningrate","text":"This function computes learning rate. It has the following inputs: outcomes: vector of 1's and 0's for behavioral outcomes windows: number of windows to split the outcome data into And the following outputs: rate: the learning rate across each window\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.mar2csd-Tuple{Any, Any, Any}","page":"API","title":"Neuroblox.mar2csd","text":"This function converts multivariate auto-regression (MAR) model parameters to a cross-spectral density (CSD). A : coefficients of MAR model, array of length p, each element contains the regression coefficients for that particular time-lag. Σ : noise covariance matrix of MAR p : number of time lags freqs : frequencies at which to evaluate the CSD sf : sampling frequency\n\nThis function returns: csd : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared)\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.mar_ml-Tuple{Any, Any}","page":"API","title":"Neuroblox.mar_ml","text":"Maximum likelihood estimator of a multivariate, or vector auto-regressive model. y : MxN Data matrix where M is number of samples and N is number of dimensions p : time lag parameter, also called order of MAR model return values mar[\"A\"] : model parameters is a NxNxP tensor, i.e. one NxN parameter matrix for each time bin k ∈ {1,...,p} mar[\"Σ\"] : noise covariance matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.matlab_norm-Tuple{Any, Any}","page":"API","title":"Neuroblox.matlab_norm","text":"function matlab_norm(A, p)\n\nSimple helper function to implement the norm of a matrix that is equivalent to the one given in MATLAB for order=1, 2, Inf. \nThis is needed for the reproduction of the exact same results of SPM12.\n\nArguments:\n- `A`: matrix\n- `p`: order of norm\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.params-Tuple{Neuroblox.BloxConnector}","page":"API","title":"Neuroblox.params","text":"Helper to merge delays and weights into a single vector\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.paramscoping-Tuple{}","page":"API","title":"Neuroblox.paramscoping","text":"function paramscoping(;kwargs...)\n\nScope arguments that are already a symbolic model parameter thereby keep the correct namespace \nand make those that are not yet symbolic a symbol.\nKeyword arguments are used, because parameter definition require names, not just values.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_cos_blox-Union{Tuple{F}, Tuple{Any, Any, F}} where F","page":"API","title":"Neuroblox.phase_cos_blox","text":"phasecosblox is creating a cos with angular frequency ω and variable phase phaseinter has the following parameters: ω: angular frequency t: time phaseinter: a function that returns phase as a function of time and returns: the resulting value\n\nUsage: phaseint = phaseinter(0:0.1:50,phasedata) phaseout(t) = phasecosblox(0.1,t,phaseint) which is now a function of time and can be used in an input blox you can also use the dot operator to calculate time-series signal = phaseout.(collect(0:0.01:50))\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_inter-Tuple{Any, Any}","page":"API","title":"Neuroblox.phase_inter","text":"phaseinter is creating a function that interpolates the phase data for any time given phaseinter has the following parameters: phaserange: a range, e.g. 0:0.1:50 which should reflect the time points of the data phasedata: phase at equidistant time points and returns: an function that returns an interpolated phase for t in range\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_sin_blox-Union{Tuple{F}, Tuple{Any, Any, F}} where F","page":"API","title":"Neuroblox.phase_sin_blox","text":"phasesinblox is creating a sin with angular frequency ω and variable phase phaseinter has the following parameters: ω: angular frequency t: time phaseinter: a function that returns phase as a function of time and returns: the resulting value\n\nUsage: phaseint = phaseinter(0:0.1:50,phasedata) phaseout(t) = phasesinblox(0.1,t,phaseint) which is now a function of time and can be used in an input blox you can also use the dot operator to calculate time-series signal = phaseout.(collect(0:0.01:50))\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phaseangle-Tuple{}","page":"API","title":"Neuroblox.phaseangle","text":"phaseangle takes in time series data, hilbert transforms it, and estimates the phase angle.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.random_initials-Tuple{ODESystem, Any}","page":"API","title":"Neuroblox.random_initials","text":"random_initials creates a vector of random initial conditions for an ODESystem that is composed of a list of blox. The function finds the initial conditions in the blox and then sets a random value in between range tuple given for that state.\n\nIt has the following inputs: odesys: ODESystem blox : list of blox\n\nAnd outputs: u0 : Float64 vector of initial conditions\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.sample_affect!-NTuple{4, Any}","page":"API","title":"Neuroblox.sample_affect!","text":"Non-symbolic, time-block-based way of `@register_symbolic sample_poisson(λ)`.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.setup_sDCM-NTuple{7, Any}","page":"API","title":"Neuroblox.setup_sDCM","text":"function setup_sDCM(data, stateevolutionmodel, initcond, csdsetup, priors, hyperpriors, indices)\n\nInterface function to performs variational inference to fit model parameters to empirical cross spectral density.\nThe current implementation provides a Variational Laplace fit (see function above `variationalbayes`).\n\nArguments:\n- `data` : dataframe with column names corresponding to the regions of measurement.\n- `model` : MTK model, including state evolution and measurement.\n- `initcond` : dictionary of initial conditions, numerical values for all states\n- `csdsetup` : dictionary of parameters required for the computation of the cross spectral density\n-- `dt` : sampling interval\n-- `freq` : frequencies at which to evaluate the CSD\n-- `p` : order parameter of the multivariate autoregression model\n- `priors` : dataframe of parameters with the following columns:\n-- `name` : corresponds to MTK model name\n-- `mean` : corresponds to prior mean value\n-- `variance` : corresponds to the prior variances\n- `hyperpriors` : dataframe of parameters with the following columns:\n-- `Πλ_pr` : prior precision matrix for λ hyperparameter(s)\n-- `μλ_pr` : prior mean(s) for λ hyperparameter(s)\n- `indices` : indices to separate model parameters from other parameters. Needed for the computation of AD gradient.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.spm_logdet-Tuple{Any}","page":"API","title":"Neuroblox.spm_logdet","text":"function spm_logdet(M)\n\nSPM12 style implementation of the logarithm of the determinant of a matrix.\n\nArguments:\n- `M`: matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.vecparam-Tuple{OrderedCollections.OrderedDict}","page":"API","title":"Neuroblox.vecparam","text":"vecparam(param::OrderedDict)\n\nFunction to flatten an ordered dictionary of model parameters and return a simple list of parameter values.\n\nArguments:\n- `param`: dictionary of model parameters (may contain numbers and lists of numbers)\n\n\n\n\n\n","category":"method"},{"location":"getting_started/#neuroblox_example","page":"Getting Started","title":"Getting Started with Neuroblox","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"This tutorial will introduce you to simulating brain dynamics using Neuroblox.","category":"page"},{"location":"getting_started/#Example-1-:-Building-an-oscillating-circuit-from-two-Wilson-Cowan-Neural-Mass-Models","page":"Getting Started","title":"Example 1 : Building an oscillating circuit from two Wilson-Cowan Neural Mass Models","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"The Wilson–Cowan model describes the dynamics of interactions between populations of excitatory and inhibitory neurons. Each Wilson-Cowan Blox is described by the follwoing equations:","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"fracdEdt = frac-Etau_E + frac11 + textexp(-a_E*(c_EE*E - c_IE*I - theta_E + eta*(sumjcn))10pt\nfracdIdt = frac-Itau_I + frac11 + exp(-a_I*(c_EI*E - c_II*I - theta_I)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Our first example is to simply combine two Wilson-Cowan Blox to build an oscillatory circuit","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"using Neuroblox\nusing DifferentialEquations\nusing Graphs\nusing MetaGraphs\nusing Plots\n\n@named WC1 = WilsonCowan()\n@named WC2 = WilsonCowan()\n\ng = MetaDiGraph()\nadd_blox!.(Ref(g), [WC1, WC2])\n\nadj = [-1 6; 6 -1]\ncreate_adjacency_edges!(g, adj)\n","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"First, we create the two Wilson-Cowan Blox: WC1 and WC2. Next, we add the two Blox into a directed graph as nodes and then we are creating weighted edges between the two nodes using an adjacency matrix.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Now we are ready to build the ModelingToolkit System. Structural simplify creates the final set of equations in which all substiutions are made.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"@named sys = system_from_graph(g)\nsys = structural_simplify(sys)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"To solve the system, we first create an ODEProblem and then solve it over the tspan of (0,100) using a stiff solver. The solution is saved every 0.1ms. The unit of time in Neuroblox is 1ms.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"prob = ODEProblem(sys, [], (0.0, 100), [])\nsol = solve(prob, Rodas4(), saveat=0.1)\nplot(sol)","category":"page"},{"location":"getting_started/#Example-2-:-Building-a-Brain-Circuit-from-literature-using-Neural-Mass-Models","page":"Getting Started","title":"Example 2 : Building a Brain Circuit from literature using Neural Mass Models","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"In this example, we will construct a Parkinsons model from eight Jansen-Rit Neural Mass Models as described in Liu et al. (2020). DOI: 10.1016/j.neunet.2019.12.021. The Jansen-Rit Neural Mass model is defined by the following differential equations:","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"fracdxdt = y-frac2taux\nfracdydt = -fracxtau^2 + fracHtau frac2lambda1+textexp(-r*sumjcn) - lambda","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"using Neuroblox\nusing DifferentialEquations\nusing Graphs\nusing MetaGraphs\nusing Plots","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"The original paper units are in seconds we therefore need to multiply all parameters with a common factor","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"τ_factor = 1000\n\n@named Str = JansenRit(τ=0.0022*τ_factor, H=20/τ_factor, λ=300, r=0.3)\n@named GPE = JansenRit(τ=0.04*τ_factor, cortical=false) # all default subcortical except τ\n@named STN = JansenRit(τ=0.01*τ_factor, H=20/τ_factor, λ=500, r=0.1)\n@named GPI = JansenRit(cortical=false) # default parameters subcortical Jansen Rit blox\n@named Th = JansenRit(τ=0.002*τ_factor, H=10/τ_factor, λ=20, r=5)\n@named EI = JansenRit(τ=0.01*τ_factor, H=20/τ_factor, λ=5, r=5)\n@named PY = JansenRit(cortical=true) # default parameters cortical Jansen Rit blox\n@named II = JansenRit(τ=2.0*τ_factor, H=60/τ_factor, λ=5, r=5)\nblox = [Str, GPE, STN, GPI, Th, EI, PY, II]","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Again, we create a graph and add the Blox as nodes","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"g = MetaDiGraph()\nadd_blox!.(Ref(g), blox)\n\nparams = @parameters C_Cor=60 C_BG_Th=60 C_Cor_BG_Th=5 C_BG_Th_Cor=5","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"ModelingToolkit allows us to create parameters that can be passed into the equations symbolically.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"We add edges as specified in Table 2 of Liu et al. We only implemented a subset of the nodes and edges to describe a less complex version of the model. Edges can also be created using an adjacency matrix as in the previous example.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"add_edge!(g, 2, 1, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 2, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 2, 3, Dict(:weight => C_BG_Th))\nadd_edge!(g, 3, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 3, 7, Dict(:weight => C_Cor_BG_Th))\nadd_edge!(g, 4, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 4, 3, Dict(:weight => C_BG_Th))\nadd_edge!(g, 5, 4, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 6, 5, Dict(:weight => C_BG_Th_Cor))\nadd_edge!(g, 6, 7, Dict(:weight => 6*C_Cor))\nadd_edge!(g, 7, 6, Dict(:weight => 4.8*C_Cor))\nadd_edge!(g, 7, 8, Dict(:weight => -1.5*C_Cor))\nadd_edge!(g, 8, 7, Dict(:weight => 1.5*C_Cor))\nadd_edge!(g, 8, 8, Dict(:weight => 3.3*C_Cor))\nadd_edge!(g,1,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,1,2,:weight, C_BG_Th)\nadd_edge!(g,2,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,2,5,:weight, C_Cor_BG_Th)\nadd_edge!(g,3,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,3,2,:weight, C_BG_Th)\nadd_edge!(g,4,3,:weight, -0.5*C_BG_Th)\nadd_edge!(g,4,4,:weight, C_BG_Th_Cor)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Now we are ready to build the ModelingToolkit System and apply structural simplification to the equations.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"@named final_system = system_from_graph(g)\nfinal_system_sys = structural_simplify(final_system)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Our Jansen-Rit model allows delayed edges, and we therefore need to collect those delays (in our case all delays are zero). Then we build a Delayed Differential Equations Problem (DDEProblem).","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"sim_dur = 1000.0 # Simulate for 1 second\nprob = ODEProblem(final_system_sys,\n [],\n (0.0, sim_dur))","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"We select an algorithm and solve the system","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"alg = Tsit5()\nsol_dde_no_delays = solve(prob, alg, saveat=1)\nplot(sol_dde_no_delays)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"In a later tutorial, we will show how to introduce edge delays.","category":"page"},{"location":"#Neuroblox","page":"Neuroblox","title":"Neuroblox","text":"","category":"section"},{"location":"#About","page":"Neuroblox","title":"About","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox.jl is designed for computational neuroscience and psychiatry applications. Our tools range from control circuit system identification to brain circuit simulations bridging scales from spiking neurons to fMRI-derived circuits, parameter-fitting models to neuroimaging data, interactions between the brain and other physiological systems, experimental optimization, and scientific machine learning.","category":"page"},{"location":"#Description","page":"Neuroblox","title":"Description","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox.jl is based on a library of modular computational building blocks (“blox”) in the form of systems of symbolic dynamic differential equations that can be combined to describe large-scale brain dynamics. Once a model is built, it can be simulated efficiently and fit electrophysiological and neuroimaging data. Moreover, the circuit behavior of multiple model variants can be investigated to aid in distinguishing between competing hypotheses. We employ ModelingToolkit.jl to describe the dynamical behavior of blox as symbolic (stochastic/delay) differential equations. Our libraries of modular blox consist of individual neurons (Hodgkin-Huxley, IF, QIF, LIF, etc.), neural mass models (Jansen-Rit, Wilson-Cowan, Lauter-Breakspear, Next Generation, microcanonical circuits etc.) and biomimetically-constrained control circuit elements. A GUI designed to be intuitive to neuroscientists allows researchers to build models that automatically generate high-performance systems of numerical ordinary/stochastic differential equations from which one can run stimulations with parameters fit to experimental data. Our benchmarks show that the increase in speed for simulation often exceeds a factor of 100 as compared to neural mass model implementation by the Virtual Brain (python) and similar packages in MATLAB. For parameter fitting of brain circuit dynamical models, we use Turing.jl to perform probabilistic modeling, including Hamilton-Monte-Carlo sampling and Automated Differentiation Variational Inference.","category":"page"},{"location":"#Installation","page":"Neuroblox","title":"Installation","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"To install Neuroblox.jl, first add the JuliaHubRegistry and then use the Julia package manager:","category":"page"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"using Pkg\nPkg.add(\"PkgAuthentication\")\nusing PkgAuthentication\nPkgAuthentication.install(\"juliahub.com\")\nPkg.Registry.add()\nPkg.add(\"Neuroblox\")","category":"page"},{"location":"#Licensing","page":"Neuroblox","title":"Licensing","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox is free for non-commerical and academic use. For full details of the license, please see the Neuroblox EULA. For commercial use, get in contact with sales@neuroblox.org.","category":"page"},{"location":"tutorials/resting_state_wb/#resting_state_tutorial","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"","category":"section"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"This tutorial will introduce you to simulating resting state brain dynamics using Neuroblox.","category":"page"},{"location":"tutorials/resting_state_wb/#Building-the-a-whole-brain-FitzHugh-Nagumo-neural-mass-model","page":"Tutorial on resting state simulation using neural mass models","title":"Building the a whole brain FitzHugh-Nagumo neural mass model","text":"","category":"section"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"The FitzHugh-Nagumo model is described by the follwoing equations:","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":" beginalign\n dotV = d tau (-f V^3 + e V^2 + alpha W - gamma I_c + sigma w(t) ) \n dotW = dfracdtau(b V - beta W + a + sigma w(t) )\n endalign","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"We start by building the resting state circuit from individual Generic2dOscillator Blox","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"using Neuroblox\nusing CSV\nusing DataFrames\nusing MetaGraphs\nusing DifferentialEquations\nusing Random\nusing Plots\nusing Statistics\nusing HypothesisTests\n\n# read connection matrix from file\nweights = CSV.read(\"../data/weights.csv\",DataFrame)\nregion_names = names(weights)\n\nwm = Array(weights)\n\n# assemble list of neural mass models\nblocks = []\nfor i in 1:size(wm)[1]\n push!(blocks, Neuroblox.Generic2dOscillator(name=Symbol(region_names[i]),bn=sqrt(5e-4)))\nend\n\n# add neural mass models to Graph and connect using the connection matrix\ng = MetaDiGraph()\nadd_blox!.(Ref(g), blocks)\ncreate_adjacency_edges!(g, wm)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"@named sys = system_from_graph(g)\nsys = structural_simplify(sys)\nlength(unknowns(sys))","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"To solve the system, we first create an Stochastic Differential Equation Problem and then solve it over the tspan of (0,6e) using a EulerHeun solver. The solution is saved every 0.5ms. The unit of time in Neuroblox is 1ms.","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"prob = SDEProblem(sys,rand(-2:0.1:4,76*2), (0.0, 6e5), [])\nsol = solve(prob, EulerHeun(), dt=0.5, saveat=5)\nplot(sol.t,sol[5,:],xlims=(0,10000))","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"To evaluate the connectivity of our simulated resting state network, we calculate the statistically significant correlations","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"cs = []\nfor i in 1:Int((length(sol.t)-1)/1000)-1\n solv = Array(sol[1:2:end,(i-1)*1000+1:(i*1000)])'\n push!(cs,cor(solv))\nend\ncss = stack(cs)\n\np = zeros(76,76)\nfor i in 1:76\n for j in 1:76\n p[i,j] = pvalue(OneSampleTTest(css[i,j,:]))\n end\nend\nheatmap(log10.(p) .* (p .< 0.05),aspect_ratio = :equal)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"Fig.: log10(p value) displaying statistally significant correlation between time series","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"heatmap(wm,aspect_ratio = :equal)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"Fig.: Connection Adjacency Matrix that was used to connect the neural mass models","category":"page"}] +[{"location":"release_notes/#Release-Notes","page":"Release Notes","title":"Release Notes","text":"","category":"section"},{"location":"release_notes/#v0.3","page":"Release Notes","title":"v0.3","text":"","category":"section"},{"location":"release_notes/","page":"Release Notes","title":"Release Notes","text":"Initial Release!","category":"page"},{"location":"api/#API-Documentation","page":"API","title":"API Documentation","text":"","category":"section"},{"location":"api/","page":"API","title":"API","text":"Modules = [Neuroblox]","category":"page"},{"location":"api/#Neuroblox.BalloonModel","page":"API","title":"Neuroblox.BalloonModel","text":"Arguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nlnκ: logarithmic prefactor to signal decay H[1], set to 0 for standard parameter value.\nlnτ: logarithmic prefactor to transit time H[3], set to 0 for standard parameter value.\nlnϵ: logarithm of ratio of intra- to extra-vascular signal\n\nNB: the prefix ln of the variables u, ν, q as well as the parameters κ, τ denotes their transformation into logarithmic space to enforce their positivity. This transformation is considered in the derivates of the model equations below. \n\nCitations:\n\nStephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040\nHofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.Generic2dOscillator","page":"API","title":"Neuroblox.Generic2dOscillator","text":"Generic2dOscillator(name, namespace, ...)\n\nThe Generic2dOscillator model is a generic dynamic system with two state\nvariables. The dynamic equations of this model are composed of two ordinary\ndifferential equations comprising two nullclines. The first nullcline is a\ncubic function as it is found in most neuron and population models; the\nsecond nullcline is arbitrarily configurable as a polynomial function up to\nsecond order. The manipulation of the latter nullcline's parameters allows\nto generate a wide range of different behaviours.\n\nEquations:\n\n```math\n \\begin{align}\n \\dot{V} &= d \\, \\tau (-f V^3 + e V^2 + g V + \\alpha W + \\gamma I) \\\\\n \\dot{W} &= \\dfrac{d}{\tau}\\,\\,(c V^2 + b V - \\beta W + a)\n \\end{align}\n```\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOther parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.\n\nCitations: FitzHugh, R., Impulses and physiological states in theoretical models of nerve membrane, Biophysical Journal 1: 445, 1961.\n\nNagumo et.al, An Active Pulse Transmission Line Simulating Nerve Axon, Proceedings of the IRE 50: 2061, 1962.\n\nStefanescu, R., Jirsa, V.K. Reduced representations of heterogeneous mixed neural networks with synaptic coupling. Physical Review E, 83, 2011.\n\nJirsa VK, Stefanescu R. Neural population modes capture biologically realistic large-scale network dynamics. Bulletin of Mathematical Biology, 2010.\n\nStefanescu, R., Jirsa, V.K. A low dimensional description of globally coupled heterogeneous neural networks of excitatory and inhibitory neurons. PLoS Computational Biology, 4(11), 2008).\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.HarmonicOscillator","page":"API","title":"Neuroblox.HarmonicOscillator","text":"HarmonicOscillator(name, namespace, ω, ζ, k, h)\n\nCreate a harmonic oscillator blox with the specified parameters.\nThe formal definition of this blox is:\n\nfracdxdt = y-(2*omega*zeta*x)+ k*(2pi)*(atan((sumjcn)h)\nfracdydt = -(omega^2)*x\n\nwhere ``jcn`` is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nω: Base frequency. Note the default value is scaled to give oscillations in milliseconds to match other blocks.\nζ: Damping ratio.\nk: Gain.\nh: Threshold.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.JansenRit","page":"API","title":"Neuroblox.JansenRit","text":"JansenRit(name, namespace, τ, H, λ, r, cortical, delayed)\n\nCreate a Jansen Rit blox as described in Liu et al.\nThe formal definition of this blox is:\n\nfracdxdt = y-frac2taux\nfracdydt = -fracxtau^2 + fracHtau frac2lambda1+textexp(-r*sumjcn) - lambda\n\nwhere jcn is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nτ: Time constant. This is changed from the original source as the time constant was in seconds, while all our blocks are in milliseconds.\nH: See equation for use.\nλ: See equation for use.\nr: See equation for use.\ncortical: Boolean to determine whether to use cortical or subcortical parameters. Specifying any of the parameters above will override this.\ndelayed: Boolean to indicate whether states are delayed\n\nCitations:\n\nLiu C, Zhou C, Wang J, Fietkiewicz C, Loparo KA. The role of coupling connections in a model of the cortico-basal ganglia-thalamocortical neural loop for the generation of beta oscillations. Neural Netw. 2020 Mar;123:381-392. doi: 10.1016/j.neunet.2019.12.021.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.JansenRitSPM12","page":"API","title":"Neuroblox.JansenRitSPM12","text":"Jansen-Rit model block for canonical micro circuit, analogous to the implementation in SPM12\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.LarterBreakspear","page":"API","title":"Neuroblox.LarterBreakspear","text":"LarterBreakspear(name, namespace, ...)\n\nCreate a Larter Breakspear blox described in Endo et al. For a full list of the parameters used see the reference.\nIf you need to modify the parameters, see Chesebro et al. and van Nieuwenhuizen et al. for physiological ranges.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOther parameters: See reference for full list. Note that parameters are scaled so that units of time are in milliseconds.\n\nCitations:\n\nEndo H, Hiroe N, Yamashita O. Evaluation of Resting Spatio-Temporal Dynamics of a Neural Mass Model Using Resting fMRI Connectivity and EEG Microstates. Front Comput Neurosci. 2020 Jan 17;13:91. doi: 10.3389/fncom.2019.00091.\nChesebro AG, Mujica-Parodi LR, Weistuch C. Ion gradient-driven bifurcations of a multi-scale neuronal model. Chaos Solitons Fractals. 2023 Feb;167:113120. doi: 10.1016/j.chaos.2023.113120. \nvan Nieuwenhuizen, H, Chesebro, AG, Polis, C, Clarke, K, Strey, HH, Weistuch, C, Mujica-Parodi, LR. Ketosis regulates K+ ion channels, strengthening brain-wide signaling disrupted by age. Preprint. bioRxiv 2023.05.10.540257; doi: https://doi.org/10.1101/2023.05.10.540257. \n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.OUBlox","page":"API","title":"Neuroblox.OUBlox","text":"Ornstein-Uhlenbeck process Blox\n\nvariables: x(t): value jcn: input parameters: τ: relaxation time \tμ: average value \tσ: random noise (variance of OU process is τ*σ^2/2) returns: an ODE System (but with brownian parameters)\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.Striatum","page":"API","title":"Neuroblox.Striatum","text":"Subcortical blox\nall subcprtical blox used in cortico-striatal model are defined here\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.WilsonCowan","page":"API","title":"Neuroblox.WilsonCowan","text":"WilsonCowan(name, namespace, τ_E, τ_I, a_E, a_I, c_EE, c_IE, c_EI, c_II, θ_E, θ_I, η)\n\nCreate a standard Wilson Cowan blox.\nThe formal definition of this blox is:\n\nfracdEdt = frac-Etau_E + frac11 + textexp(-a_E*(c_EE*E - c_IE*I - theta_E + eta*(sumjcn))\nfracdIdt = frac-Itau_I + frac11 + exp(-a_I*(c_EI*E - c_II*I - theta_I)\n\nwhere jcn is any input to the blox.\n\nArguments:\n\nname: Name given to ODESystem object within the blox.\nnamespace: Additional namespace above name if needed for inheritance.\nOthers: See equation for use.\n\n\n\n\n\n","category":"type"},{"location":"api/#Neuroblox.WinnerTakeAllBlox","page":"API","title":"Neuroblox.WinnerTakeAllBlox","text":"WinnerTakeAllBlox\n\nCreates a winner-take-all local circuit found in neocortex, typically 5 pyramidal (excitatory) neurons send synapses to a single interneuron (inhibitory) and receive feedback inhibition from that interneuron.\n\n\n\n\n\n","category":"type"},{"location":"api/#LinearAlgebra.eigen-Union{Tuple{Array{ForwardDiff.Dual{T, P, np}, 2}}, Tuple{np}, Tuple{P}, Tuple{T}} where {T, P, np}","page":"API","title":"LinearAlgebra.eigen","text":"function LinearAlgebra.eigen(M::Matrix{Dual{T, P, np}}) where {T, P, np}\n\nDispatch of LinearAlgebra.eigen for dual matrices with complex numbers. Make the eigenvalue decomposition \namenable to automatic differentiation. To do so compute the analytical derivative of eigenvalues\nand eigenvectors. \n\nArguments:\n- `M`: matrix of type Dual of which to compute the eigenvalue decomposition. \n\nReturns:\n- `Eigen(evals, evecs)`: eigenvalue decomposition returned as type LinearAlgebra.Eigen\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.ARVTarget-NTuple{6, Any}","page":"API","title":"Neuroblox.ARVTarget","text":"ARVTarget Time series data is bandpass filtered and then the power spectrum is computed for a given time interval (control bin), returned as the average value of the power spectral density within a certain frequency band ([lb, ub]).\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.CDVTarget-NTuple{5, Any}","page":"API","title":"Neuroblox.CDVTarget","text":"CDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Circular difference is quantified as the angle of circular_location.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.ControlError-NTuple{8, Any}","page":"API","title":"Neuroblox.ControlError","text":"ControlError Returns the control error (deviation of the actual value from the target value).\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.PDVTarget-NTuple{5, Any}","page":"API","title":"Neuroblox.PDVTarget","text":"PDVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians. Phase deviation is quantified as the angle difference between a given set of signals.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.PLVTarget-NTuple{6, Any}","page":"API","title":"Neuroblox.PLVTarget","text":"PLVTarget Time series data is bandpass filtered and hilbert-transformed. Phase angle is computed in radians.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.addnontunableparams-Tuple{Any, Any}","page":"API","title":"Neuroblox.addnontunableparams","text":"function addnontunableparams(param, model)\n\nFunction adds parameters of a model that were not marked as tunable to a list of tunable parameters\nand respects the MTK ordering of parameters.\n\nArguments:\n- `paramlist`: parameters of an MTK system that were tagged as tunable\n- `sys`: MTK system\n\nReturns:\n- `completeparamlist`: complete parameter list of a system, including those that were not tagged as tunable\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.bandpassfilter-Tuple{}","page":"API","title":"Neuroblox.bandpassfilter","text":"bandpassfilter takes in time series data and bandpass filters it. It has the following inputs: data: time series data lb: minimum cut-off frequency ub: maximum cut-off frequency fs: sampling frequency order: filter order\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.boldsignal-Tuple{}","page":"API","title":"Neuroblox.boldsignal","text":"Arguments:\n\nname: Name given to ODESystem object within the blox.\nlnϵ : logarithm of ratio of intra- to extra-vascular signal\n\nNB: the prefix ln of the variables ν, q as well as the parameters ϵ denotes their transformation into logarithmic space to enforce their positivity.\n\nCitations:\n\nStephan K E, Weiskopf N, Drysdale P M, Robinson P A, and Friston K J. Comparing Hemodynamic Models with DCM. NeuroImage 38, no. 3 (2007): 387–401. doi: 10.1016/j.neuroimage.2007.07.040\nHofmann D, Chesebro A G, Rackauckas C, Mujica-Parodi L R, Friston K J, Edelman A, and Strey H H. Leveraging Julia's Automated Differentiation and Symbolic Computation to Increase Spectral DCM Flexibility and Speed, 2023. doi: 10.1101/2023.10.27.564407\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.complexwavelet","page":"API","title":"Neuroblox.complexwavelet","text":"complexwavelet creates a complex morlet wavelet by windowing a complex sine wave with a Gaussian taper. The morlet wavelet is a special case of a bandpass filter in which the frequency response is Gaussian-shaped. Convolution with a complex wavelet is equivalent to performing a Hilbert transform of a bandpass filtered signal.\n\nIt has the following inputs: data: time series data dt : data sampling rate lb : lower bound wavelet frequency (in Hz) ub : upper bound wavelet frequency (in Hz) a : amplitude of the Gaussian taper, default is 1 n : number of wavelet cycles of the Gaussian taper, defines the trade-off between temporal precision and frequency precision larger n gives better frequency precision at the cost of temporal precision default is 6 Hz m : x-axis offset, default is 0 num_wavelets : number of wavelets to create, default is 5\n\nAnd outputs: complex_wavelet : a family of complex morlet wavelets\n\n\n\n\n\n","category":"function"},{"location":"api/#Neuroblox.csd2mar-NTuple{4, Any}","page":"API","title":"Neuroblox.csd2mar","text":"This function converts a cross-spectral density (CSD) into a multivariate auto-regression (MAR) model. It first transforms the CSD into its cross-correlation function (Wiener-Kinchine theorem) and then computes the MAR model coefficients. csd : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared) w : frequencies dt : time step size p : number of time steps of auto-regressive model\n\nThis function returns coeff : array of length p of coefficient matrices of size sqrt(N)xsqrt(N) noise_cov : noise covariance matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.csd_approx-NTuple{4, Any}","page":"API","title":"Neuroblox.csd_approx","text":"This function implements equation 2 of the spectral DCM paper, Friston et al. 2014 \"A DCM for resting state fMRI\".\nNote that nomenclature is taken from SPM12 code and it does not seem to coincide with the spectral DCM paper's nomenclature. \nFor instance, Gu should represent the spectral component due to external input according to the paper. However, in the code this represents\nthe hidden state fluctuations (which are called Gν in the paper).\nGn in the code corresponds to Ge in the paper, i.e. the observation noise. In the code global and local components are defined, no such distinction\nis discussed in the paper. In fact the parameter γ, corresponding to local component is not present in the paper.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.get_dynamic_states-Tuple{Any}","page":"API","title":"Neuroblox.get_dynamic_states","text":"function get_dynamic_states(sys)\n\nFunction extracts states from the system that are dynamic variables, \nget also indices of external inputs (u(t)) and measurements (like bold(t))\nArguments:\n- `sys`: MTK system\n\nReturns:\n- `sts`: states/unknowns of the system that are neither external inputs nor measurements, i.e. these are the dynamic states\n- `idx`: indices of these states\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.idft-Tuple{AbstractArray}","page":"API","title":"Neuroblox.idft","text":"Plain implementation of idft because AD dispatch versions for ifft don't work still!\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.inner_namespaceof-Tuple{Any}","page":"API","title":"Neuroblox.inner_namespaceof","text":"Returns the complete namespace EXCLUDING the outermost (highest) level.\nThis is useful for manually preparing equations (e.g. connections, see BloxConnector),\nthat will later be composed and will automatically get the outermost namespace.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.input_equations-Tuple{Any}","page":"API","title":"Neuroblox.input_equations","text":"Returns the equations for all input variables of a system, \nassuming they have a form like : `sys.input_variable ~ ...`\nso only the input appears on the LHS.\n\nInput equations are namespaced by the inner namespace of blox\nand then they are returned. This way during system `compose` downstream,\nthe higher-level namespaces will be added to them.\n\nIf blox isa AbstractComponent, it is assumed that it contains a `connector` field,\nwhich holds a `BloxConnector` object with all relevant connections \nfrom lower levels and this level.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.learningrate-Tuple{Any, Any}","page":"API","title":"Neuroblox.learningrate","text":"This function computes learning rate. It has the following inputs: outcomes: vector of 1's and 0's for behavioral outcomes windows: number of windows to split the outcome data into And the following outputs: rate: the learning rate across each window\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.mar2csd-Tuple{Any, Any, Any}","page":"API","title":"Neuroblox.mar2csd","text":"This function converts multivariate auto-regression (MAR) model parameters to a cross-spectral density (CSD). A : coefficients of MAR model, array of length p, each element contains the regression coefficients for that particular time-lag. Σ : noise covariance matrix of MAR p : number of time lags freqs : frequencies at which to evaluate the CSD sf : sampling frequency\n\nThis function returns: csd : cross-spectral density matrix of size MxN; M: number of samples, N: number of cross-spectral dimensions (number of variables squared)\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.mar_ml-Tuple{Any, Any}","page":"API","title":"Neuroblox.mar_ml","text":"Maximum likelihood estimator of a multivariate, or vector auto-regressive model. y : MxN Data matrix where M is number of samples and N is number of dimensions p : time lag parameter, also called order of MAR model return values mar[\"A\"] : model parameters is a NxNxP tensor, i.e. one NxN parameter matrix for each time bin k ∈ {1,...,p} mar[\"Σ\"] : noise covariance matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.matlab_norm-Tuple{Any, Any}","page":"API","title":"Neuroblox.matlab_norm","text":"function matlab_norm(A, p)\n\nSimple helper function to implement the norm of a matrix that is equivalent to the one given in MATLAB for order=1, 2, Inf. \nThis is needed for the reproduction of the exact same results of SPM12.\n\nArguments:\n- `A`: matrix\n- `p`: order of norm\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.params-Tuple{Neuroblox.BloxConnector}","page":"API","title":"Neuroblox.params","text":"Helper to merge delays and weights into a single vector\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.paramscoping-Tuple{}","page":"API","title":"Neuroblox.paramscoping","text":"function paramscoping(;kwargs...)\n\nScope arguments that are already a symbolic model parameter thereby keep the correct namespace \nand make those that are not yet symbolic a symbol.\nKeyword arguments are used, because parameter definition require names, not just values.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_cos_blox-Union{Tuple{F}, Tuple{Any, Any, F}} where F","page":"API","title":"Neuroblox.phase_cos_blox","text":"phasecosblox is creating a cos with angular frequency ω and variable phase phaseinter has the following parameters: ω: angular frequency t: time phaseinter: a function that returns phase as a function of time and returns: the resulting value\n\nUsage: phaseint = phaseinter(0:0.1:50,phasedata) phaseout(t) = phasecosblox(0.1,t,phaseint) which is now a function of time and can be used in an input blox you can also use the dot operator to calculate time-series signal = phaseout.(collect(0:0.01:50))\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_inter-Tuple{Any, Any}","page":"API","title":"Neuroblox.phase_inter","text":"phaseinter is creating a function that interpolates the phase data for any time given phaseinter has the following parameters: phaserange: a range, e.g. 0:0.1:50 which should reflect the time points of the data phasedata: phase at equidistant time points and returns: an function that returns an interpolated phase for t in range\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phase_sin_blox-Union{Tuple{F}, Tuple{Any, Any, F}} where F","page":"API","title":"Neuroblox.phase_sin_blox","text":"phasesinblox is creating a sin with angular frequency ω and variable phase phaseinter has the following parameters: ω: angular frequency t: time phaseinter: a function that returns phase as a function of time and returns: the resulting value\n\nUsage: phaseint = phaseinter(0:0.1:50,phasedata) phaseout(t) = phasesinblox(0.1,t,phaseint) which is now a function of time and can be used in an input blox you can also use the dot operator to calculate time-series signal = phaseout.(collect(0:0.01:50))\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.phaseangle-Tuple{}","page":"API","title":"Neuroblox.phaseangle","text":"phaseangle takes in time series data, hilbert transforms it, and estimates the phase angle.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.random_initials-Tuple{ODESystem, Any}","page":"API","title":"Neuroblox.random_initials","text":"random_initials creates a vector of random initial conditions for an ODESystem that is composed of a list of blox. The function finds the initial conditions in the blox and then sets a random value in between range tuple given for that state.\n\nIt has the following inputs: odesys: ODESystem blox : list of blox\n\nAnd outputs: u0 : Float64 vector of initial conditions\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.sample_affect!-NTuple{4, Any}","page":"API","title":"Neuroblox.sample_affect!","text":"Non-symbolic, time-block-based way of `@register_symbolic sample_poisson(λ)`.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.setup_sDCM-NTuple{7, Any}","page":"API","title":"Neuroblox.setup_sDCM","text":"function setup_sDCM(data, stateevolutionmodel, initcond, csdsetup, priors, hyperpriors, indices)\n\nInterface function to performs variational inference to fit model parameters to empirical cross spectral density.\nThe current implementation provides a Variational Laplace fit (see function above `variationalbayes`).\n\nArguments:\n- `data` : dataframe with column names corresponding to the regions of measurement.\n- `model` : MTK model, including state evolution and measurement.\n- `initcond` : dictionary of initial conditions, numerical values for all states\n- `csdsetup` : dictionary of parameters required for the computation of the cross spectral density\n-- `dt` : sampling interval\n-- `freq` : frequencies at which to evaluate the CSD\n-- `p` : order parameter of the multivariate autoregression model\n- `priors` : dataframe of parameters with the following columns:\n-- `name` : corresponds to MTK model name\n-- `mean` : corresponds to prior mean value\n-- `variance` : corresponds to the prior variances\n- `hyperpriors` : dataframe of parameters with the following columns:\n-- `Πλ_pr` : prior precision matrix for λ hyperparameter(s)\n-- `μλ_pr` : prior mean(s) for λ hyperparameter(s)\n- `indices` : indices to separate model parameters from other parameters. Needed for the computation of AD gradient.\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.spm_logdet-Tuple{Any}","page":"API","title":"Neuroblox.spm_logdet","text":"function spm_logdet(M)\n\nSPM12 style implementation of the logarithm of the determinant of a matrix.\n\nArguments:\n- `M`: matrix\n\n\n\n\n\n","category":"method"},{"location":"api/#Neuroblox.vecparam-Tuple{OrderedCollections.OrderedDict}","page":"API","title":"Neuroblox.vecparam","text":"vecparam(param::OrderedDict)\n\nFunction to flatten an ordered dictionary of model parameters and return a simple list of parameter values.\n\nArguments:\n- `param`: dictionary of model parameters (may contain numbers and lists of numbers)\n\n\n\n\n\n","category":"method"},{"location":"getting_started/#neuroblox_example","page":"Getting Started","title":"Getting Started with Neuroblox","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"This tutorial will introduce you to simulating brain dynamics using Neuroblox.","category":"page"},{"location":"getting_started/#Example-1-:-Building-an-oscillating-circuit-from-two-Wilson-Cowan-Neural-Mass-Models","page":"Getting Started","title":"Example 1 : Building an oscillating circuit from two Wilson-Cowan Neural Mass Models","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"The Wilson–Cowan model describes the dynamics of interactions between populations of excitatory and inhibitory neurons. Each Wilson-Cowan Blox is described by the follwoing equations:","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"fracdEdt = frac-Etau_E + frac11 + textexp(-a_E*(c_EE*E - c_IE*I - theta_E + eta*(sumjcn))10pt\nfracdIdt = frac-Itau_I + frac11 + exp(-a_I*(c_EI*E - c_II*I - theta_I)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Our first example is to simply combine two Wilson-Cowan Blox to build an oscillatory circuit","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"using Neuroblox\nusing DifferentialEquations\nusing Graphs\nusing MetaGraphs\nusing Plots\n\n@named WC1 = WilsonCowan()\n@named WC2 = WilsonCowan()\n\ng = MetaDiGraph()\nadd_blox!.(Ref(g), [WC1, WC2])\n\nadj = [-1 6; 6 -1]\ncreate_adjacency_edges!(g, adj)\n","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"First, we create the two Wilson-Cowan Blox: WC1 and WC2. Next, we add the two Blox into a directed graph as nodes and then we are creating weighted edges between the two nodes using an adjacency matrix.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Now we are ready to build the ModelingToolkit System. Structural simplify creates the final set of equations in which all substiutions are made.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"@named sys = system_from_graph(g)\nsys = structural_simplify(sys)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"To solve the system, we first create an ODEProblem and then solve it over the tspan of (0,100) using a stiff solver. The solution is saved every 0.1ms. The unit of time in Neuroblox is 1ms.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"prob = ODEProblem(sys, [], (0.0, 100), [])\nsol = solve(prob, Rodas4(), saveat=0.1)\nplot(sol)","category":"page"},{"location":"getting_started/#Example-2-:-Building-a-Brain-Circuit-from-literature-using-Neural-Mass-Models","page":"Getting Started","title":"Example 2 : Building a Brain Circuit from literature using Neural Mass Models","text":"","category":"section"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"In this example, we will construct a Parkinsons model from eight Jansen-Rit Neural Mass Models as described in Liu et al. (2020). DOI: 10.1016/j.neunet.2019.12.021. The Jansen-Rit Neural Mass model is defined by the following differential equations:","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"fracdxdt = y-frac2taux\nfracdydt = -fracxtau^2 + fracHtau frac2lambda1+textexp(-r*sumjcn) - lambda","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"using Neuroblox\nusing DifferentialEquations\nusing Graphs\nusing MetaGraphs\nusing Plots","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"The original paper units are in seconds we therefore need to multiply all parameters with a common factor","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"τ_factor = 1000\n\n@named Str = JansenRit(τ=0.0022*τ_factor, H=20/τ_factor, λ=300, r=0.3)\n@named GPE = JansenRit(τ=0.04*τ_factor, cortical=false) # all default subcortical except τ\n@named STN = JansenRit(τ=0.01*τ_factor, H=20/τ_factor, λ=500, r=0.1)\n@named GPI = JansenRit(cortical=false) # default parameters subcortical Jansen Rit blox\n@named Th = JansenRit(τ=0.002*τ_factor, H=10/τ_factor, λ=20, r=5)\n@named EI = JansenRit(τ=0.01*τ_factor, H=20/τ_factor, λ=5, r=5)\n@named PY = JansenRit(cortical=true) # default parameters cortical Jansen Rit blox\n@named II = JansenRit(τ=2.0*τ_factor, H=60/τ_factor, λ=5, r=5)\nblox = [Str, GPE, STN, GPI, Th, EI, PY, II]","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Again, we create a graph and add the Blox as nodes","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"g = MetaDiGraph()\nadd_blox!.(Ref(g), blox)\n\nparams = @parameters C_Cor=60 C_BG_Th=60 C_Cor_BG_Th=5 C_BG_Th_Cor=5","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"ModelingToolkit allows us to create parameters that can be passed into the equations symbolically.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"We add edges as specified in Table 2 of Liu et al. We only implemented a subset of the nodes and edges to describe a less complex version of the model. Edges can also be created using an adjacency matrix as in the previous example.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"add_edge!(g, 2, 1, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 2, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 2, 3, Dict(:weight => C_BG_Th))\nadd_edge!(g, 3, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 3, 7, Dict(:weight => C_Cor_BG_Th))\nadd_edge!(g, 4, 2, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 4, 3, Dict(:weight => C_BG_Th))\nadd_edge!(g, 5, 4, Dict(:weight => -0.5*C_BG_Th))\nadd_edge!(g, 6, 5, Dict(:weight => C_BG_Th_Cor))\nadd_edge!(g, 6, 7, Dict(:weight => 6*C_Cor))\nadd_edge!(g, 7, 6, Dict(:weight => 4.8*C_Cor))\nadd_edge!(g, 7, 8, Dict(:weight => -1.5*C_Cor))\nadd_edge!(g, 8, 7, Dict(:weight => 1.5*C_Cor))\nadd_edge!(g, 8, 8, Dict(:weight => 3.3*C_Cor))\nadd_edge!(g,1,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,1,2,:weight, C_BG_Th)\nadd_edge!(g,2,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,2,5,:weight, C_Cor_BG_Th)\nadd_edge!(g,3,1,:weight, -0.5*C_BG_Th)\nadd_edge!(g,3,2,:weight, C_BG_Th)\nadd_edge!(g,4,3,:weight, -0.5*C_BG_Th)\nadd_edge!(g,4,4,:weight, C_BG_Th_Cor)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Now we are ready to build the ModelingToolkit System and apply structural simplification to the equations.","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"@named final_system = system_from_graph(g)\nfinal_system_sys = structural_simplify(final_system)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"Our Jansen-Rit model allows delayed edges, and we therefore need to collect those delays (in our case all delays are zero). Then we build a Delayed Differential Equations Problem (DDEProblem).","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"sim_dur = 1000.0 # Simulate for 1 second\nprob = ODEProblem(final_system_sys,\n [],\n (0.0, sim_dur))","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"We select an algorithm and solve the system","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"alg = Tsit5()\nsol_dde_no_delays = solve(prob, alg, saveat=1)\nplot(sol_dde_no_delays)","category":"page"},{"location":"getting_started/","page":"Getting Started","title":"Getting Started","text":"In a later tutorial, we will show how to introduce edge delays.","category":"page"},{"location":"#Neuroblox","page":"Neuroblox","title":"Neuroblox","text":"","category":"section"},{"location":"#About","page":"Neuroblox","title":"About","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox.jl is designed for computational neuroscience and psychiatry applications. Our tools range from control circuit system identification to brain circuit simulations bridging scales from spiking neurons to fMRI-derived circuits, parameter-fitting models to neuroimaging data, interactions between the brain and other physiological systems, experimental optimization, and scientific machine learning.","category":"page"},{"location":"#Description","page":"Neuroblox","title":"Description","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox.jl is based on a library of modular computational building blocks (“blox”) in the form of systems of symbolic dynamic differential equations that can be combined to describe large-scale brain dynamics. Once a model is built, it can be simulated efficiently and fit electrophysiological and neuroimaging data. Moreover, the circuit behavior of multiple model variants can be investigated to aid in distinguishing between competing hypotheses. We employ ModelingToolkit.jl to describe the dynamical behavior of blox as symbolic (stochastic/delay) differential equations. Our libraries of modular blox consist of individual neurons (Hodgkin-Huxley, IF, QIF, LIF, etc.), neural mass models (Jansen-Rit, Wilson-Cowan, Lauter-Breakspear, Next Generation, microcanonical circuits etc.) and biomimetically-constrained control circuit elements. A GUI designed to be intuitive to neuroscientists allows researchers to build models that automatically generate high-performance systems of numerical ordinary/stochastic differential equations from which one can run stimulations with parameters fit to experimental data. Our benchmarks show that the increase in speed for simulation often exceeds a factor of 100 as compared to neural mass model implementation by the Virtual Brain (python) and similar packages in MATLAB. For parameter fitting of brain circuit dynamical models, we use Turing.jl to perform probabilistic modeling, including Hamilton-Monte-Carlo sampling and Automated Differentiation Variational Inference.","category":"page"},{"location":"#Installation","page":"Neuroblox","title":"Installation","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"To install Neuroblox.jl, first add the JuliaHubRegistry and then use the Julia package manager:","category":"page"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"using Pkg\nPkg.add(\"PkgAuthentication\")\nusing PkgAuthentication\nPkgAuthentication.install(\"juliahub.com\")\nPkg.Registry.add()\nPkg.add(\"Neuroblox\")","category":"page"},{"location":"#Licensing","page":"Neuroblox","title":"Licensing","text":"","category":"section"},{"location":"","page":"Neuroblox","title":"Neuroblox","text":"Neuroblox is free for non-commerical and academic use. For full details of the license, please see the Neuroblox EULA. For commercial use, get in contact with sales@neuroblox.org.","category":"page"},{"location":"tutorials/resting_state_wb/#resting_state_tutorial","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"","category":"section"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"This tutorial will introduce you to simulating resting state brain dynamics using Neuroblox.","category":"page"},{"location":"tutorials/resting_state_wb/#Building-the-a-whole-brain-FitzHugh-Nagumo-neural-mass-model","page":"Tutorial on resting state simulation using neural mass models","title":"Building the a whole brain FitzHugh-Nagumo neural mass model","text":"","category":"section"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"The FitzHugh-Nagumo model is described by the follwoing equations:","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":" beginalign\n dotV = d tau (-f V^3 + e V^2 + alpha W - gamma I_c + sigma w(t) ) \n dotW = dfracdtau(b V - beta W + a + sigma w(t) )\n endalign","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"We start by building the resting state circuit from individual Generic2dOscillator Blox","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"using Neuroblox\nusing CSV\nusing DataFrames\nusing MetaGraphs\nusing DifferentialEquations\nusing Random\nusing Plots\nusing Statistics\nusing HypothesisTests\n\n# read connection matrix from file\nweights = CSV.read(\"../data/weights.csv\",DataFrame)\nregion_names = names(weights)\n\nwm = Array(weights)\n\n# assemble list of neural mass models\nblocks = []\nfor i in 1:size(wm)[1]\n push!(blocks, Neuroblox.Generic2dOscillator(name=Symbol(region_names[i]),bn=sqrt(5e-4)))\nend\n\n# add neural mass models to Graph and connect using the connection matrix\ng = MetaDiGraph()\nadd_blox!.(Ref(g), blocks)\ncreate_adjacency_edges!(g, wm)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"@named sys = system_from_graph(g)\nsys = structural_simplify(sys)\nlength(unknowns(sys))","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"To solve the system, we first create an Stochastic Differential Equation Problem and then solve it over the tspan of (0,6e) using a EulerHeun solver. The solution is saved every 0.5ms. The unit of time in Neuroblox is 1ms.","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"prob = SDEProblem(sys,rand(-2:0.1:4,76*2), (0.0, 6e5), [])\nsol = solve(prob, EulerHeun(), dt=0.5, saveat=5)\nplot(sol.t,sol[5,:],xlims=(0,10000))","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"To evaluate the connectivity of our simulated resting state network, we calculate the statistically significant correlations","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"cs = []\nfor i in 1:Int((length(sol.t)-1)/1000)-1\n solv = Array(sol[1:2:end,(i-1)*1000+1:(i*1000)])'\n push!(cs,cor(solv))\nend\ncss = stack(cs)\n\np = zeros(76,76)\nfor i in 1:76\n for j in 1:76\n p[i,j] = pvalue(OneSampleTTest(css[i,j,:]))\n end\nend\nheatmap(log10.(p) .* (p .< 0.05),aspect_ratio = :equal)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"Fig.: log10(p value) displaying statistally significant correlation between time series","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"heatmap(wm,aspect_ratio = :equal)","category":"page"},{"location":"tutorials/resting_state_wb/","page":"Tutorial on resting state simulation using neural mass models","title":"Tutorial on resting state simulation using neural mass models","text":"Fig.: Connection Adjacency Matrix that was used to connect the neural mass models","category":"page"}] } diff --git a/dev/tutorials/resting_state_wb/4cba14ce.svg b/dev/tutorials/resting_state_wb/4cba14ce.svg new file mode 100644 index 0000000..4332536 --- /dev/null +++ b/dev/tutorials/resting_state_wb/4cba14ce.svg @@ -0,0 +1,46 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/dev/tutorials/resting_state_wb/78ec20c5.svg b/dev/tutorials/resting_state_wb/784ab2a6.svg similarity index 90% rename from dev/tutorials/resting_state_wb/78ec20c5.svg rename to dev/tutorials/resting_state_wb/784ab2a6.svg index f237f49..d4c19bb 100644 --- a/dev/tutorials/resting_state_wb/78ec20c5.svg +++ b/dev/tutorials/resting_state_wb/784ab2a6.svg @@ -1,49 +1,49 @@ - + - + - + - + - + - - - - - - - - - - - - - - - - - - - - - - - - - - - + + + + + + + + + + + + + + + + + + + + + + + + + + + - + - + - - - - - - - - + + + + + + + + diff --git a/dev/tutorials/resting_state_wb/b0e80753.svg b/dev/tutorials/resting_state_wb/b0e80753.svg new file mode 100644 index 0000000..d2cb9bf --- /dev/null +++ b/dev/tutorials/resting_state_wb/b0e80753.svg @@ -0,0 +1,547 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/dev/tutorials/resting_state_wb/d8b495bd.svg b/dev/tutorials/resting_state_wb/d8b495bd.svg deleted file mode 100644 index 9659d2a..0000000 --- a/dev/tutorials/resting_state_wb/d8b495bd.svg +++ /dev/null @@ -1,551 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - diff --git a/dev/tutorials/resting_state_wb/e5d28944.svg b/dev/tutorials/resting_state_wb/e5d28944.svg deleted file mode 100644 index 8d107be..0000000 --- a/dev/tutorials/resting_state_wb/e5d28944.svg +++ /dev/null @@ -1,48 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - diff --git a/dev/tutorials/resting_state_wb/index.html b/dev/tutorials/resting_state_wb/index.html index dd342fa..dece3dc 100644 --- a/dev/tutorials/resting_state_wb/index.html +++ b/dev/tutorials/resting_state_wb/index.html @@ -31,7 +31,7 @@ sys = structural_simplify(sys) length(unknowns(sys))
152

To solve the system, we first create an Stochastic Differential Equation Problem and then solve it over the tspan of (0,6e) using a EulerHeun solver. The solution is saved every 0.5ms. The unit of time in Neuroblox is 1ms.

prob = SDEProblem(sys,rand(-2:0.1:4,76*2), (0.0, 6e5), [])
 sol = solve(prob, EulerHeun(), dt=0.5, saveat=5)
-plot(sol.t,sol[5,:],xlims=(0,10000))
Example block output

To evaluate the connectivity of our simulated resting state network, we calculate the statistically significant correlations

cs = []
+plot(sol.t,sol[5,:],xlims=(0,10000))
Example block output

To evaluate the connectivity of our simulated resting state network, we calculate the statistically significant correlations

cs = []
 for i in 1:Int((length(sol.t)-1)/1000)-1
     solv = Array(sol[1:2:end,(i-1)*1000+1:(i*1000)])'
     push!(cs,cor(solv))
@@ -44,4 +44,4 @@
         p[i,j] = pvalue(OneSampleTTest(css[i,j,:]))
     end
 end
-heatmap(log10.(p) .* (p .< 0.05),aspect_ratio = :equal)
Example block output

Fig.: log10(p value) displaying statistally significant correlation between time series

heatmap(wm,aspect_ratio = :equal)
Example block output

Fig.: Connection Adjacency Matrix that was used to connect the neural mass models

+heatmap(log10.(p) .* (p .< 0.05),aspect_ratio = :equal)Example block output

Fig.: log10(p value) displaying statistally significant correlation between time series

heatmap(wm,aspect_ratio = :equal)
Example block output

Fig.: Connection Adjacency Matrix that was used to connect the neural mass models