-
Notifications
You must be signed in to change notification settings - Fork 6
Home
- The COIN class
- Running a COIN model simulation
- Performing multiple runs of a COIN model simulation
- Plotting COIN model variables
- Fitting the parameters of the COIN model to data
The COIN model is implemented as a class in MATLAB called COIN. An object of the class can be created by calling the class name:
obj = COIN
This object has properties that define the model (e.g., model parameters) and the paradigm (e.g., perturbations, sensory cues). Some properties have default values (set at the top of COIN.m), which can be overwritten. Data stored in properties can be operated on by methods (functions) of the class.
For a full list of properties of the COIN class, call help COIN (see 'PROPERTIES' section).
To run a COIN model simulation, first define a trial-by-trial sequence of perturbations via the 'perturbations' property (use NaN to indicate a channel trial). For example, for a spontaneous recovery paradigm:
obj.perturbations = [zeros(1,50) ones(1,125) -ones(1,15) NaN(1,150)];
If an experiment involves sensory cues, these can also be provided via the 'cues' property. Next call the simulate_COIN method using the dot notation:
OUTPUT = obj.simulate_COIN;
The simulation results are contained in the structure OUTPUT. The state feedback and motor output of the COIN model can be plotted as follows:
figure
plot(OUTPUT.runs{1}.state_feedback,'.')
hold on
plot(OUTPUT.runs{1}.motor_output)
legend('state feedback','motor output')
The state feedback observations are generated randomly within the simulate_COIN method by adding observation noise (sampled from the prior) to the user-defined perturbations.
Multiple calls to the simulate_COIN method will produce different results, even for the same set of parameters and perturbation sequence, as inference in the COIN model depends on state feedback, which is generated randomly within the simulate_COIN method. Such variability in model behaviour is often undesirable. To produce consistent model behaviour for a given set of parameters and perturbation sequence, multiple runs of a simulation can be performed, each conditioned on a different state feedback sequence, and the results of the runs averaged. In general, averaging over more runs will produce cleaner, less variable results.
To specify the number of runs to perform, use the 'runs' property. For example, to perform 10 runs:
obj = COIN
obj.perturbations = [zeros(1,50) ones(1,125) -ones(1,15) NaN(1,150)];
obj.runs = 10;
OUTPUT = obj.simulate_COIN;
OUTPUT.runs is a cell array with one cell per run and OUTPUT.weights is a vector specifying the relative weight of each run (weights are non-negative and sum to 1). These weights should be used to compute a weighted average of the motor output of the COIN model over runs:
for run = 1:obj.runs
motor_output(run,:) = OUTPUT.runs{run}.motor_output;
end
plot(OUTPUT.weights*motor_output)
In this example, each run is assigned an equal weight, as no adaptation data was passed to the class object. However, if the parameters of the COIN model have been first fit to adaptation data, and this data is passed to the class object before calling the simulate_COIN method, each run will be assigned a weight based on the likelihood of the adaptation data (see Performing multiple runs of a COIN model simulation conditioned on the adaptation data of a participant).
The time it takes to execute multiple runs of a COIN model simulation in series (e.g., using a for loop) can be prohibitively long if the number of runs is large. To execute multiple runs in parallel across different CPU cores and thus reduce the model runtime, specify the maximum number of CPU cores available using the 'max_cores' property (the default setting of max_cores is 0, which executes multiple runs in series). You may want to use a computer cluster.
The distributions or expected values of variables inferred by the COIN model can be plotted by activating the relevant properties of the class object. For example, to plot the predicted state distribution for each context, the predicted probabilities and the overall predicted state distribution:
obj = COIN
obj.perturbations = [zeros(1,50) ones(1,125) -ones(1,15) NaN(1,150)];
obj.runs = 10;
obj.max_cores = feature('numcores');
obj.plot_state_given_context = true;
obj.plot_predicted_probabilities = true;
obj.plot_state = true;
obj.state_values = linspace(-1.5,1.5,500); % points to evaluate the ‘state | context’ distribution at
OUTPUT = obj.simulate_COIN;
Note that for distributions of continuous variables, the points at which to evaluate the distribution can be specified (as shown in the example above). Next call the plot_COIN method:
PLOTS = obj.plot_COIN(OUTPUT);
This will generate the requested plots—a state | context plot and a predicted probabilities plot in this example. The structure PLOTS contains the data that is plotted (to see how the data is plotted, view the generate_figures method in COIN.m). Note that these plots may take some time to generate, as they require contexts to be relabelled in multiple particles and multiple runs on each trial. Once the contexts have been relabelled, the relevant distributions or expected values of variables are automatically averaged over particles and runs. In general, using more runs will result in cleaner, less variable plots. For a full list of variables that can be plotted and how to plot them, call help COIN, and see the 'plot flags' and 'plot inputs' subsections of the 'PROPERTIES' section.
The parameters of the COIN model are fit to data by minimising the negative log of the likelihood function. For an optimiser to perform this minimisation, it must have access to a function that can evaluate the negative log-likelihood as a function of the model parameters. In the COIN model, this function is performed by the objective_COIN method:
negative_log_likelihood = obj.objective_COIN;
The objective_COIN method returns a stochastic estimate of the negative log-likelihood. The estimate is stochastic because it is calculated by simulating one or more runs of the COIN model, each conditioned on a different random sequence of state feedback observations. To aid parameter optimisation, the variance of the estimate of the negative log-likelihood can be reduced by increasing the number of runs used to calculate the negative log-likelihood via the 'runs' property; to avoid excessive runtimes, multiple runs should be performed in parallel. Because the estimate of the negative log-likelihood is stochastic, an optimiser that can handle a stochastic objective function should be used to fit the parameters of the COIN model to data (e.g., BADS).
To view code showing how to fit the parameters of the COIN model to the data of an individual participant, see fit_parameters_to_data.m (BADS installation required).
Performing multiple runs of a COIN model simulation conditioned on the adaptation data of a participant
Once the parameters of the COIN model have been fit to the adaptation data of a participant, each run of a COIN model simulation can be assigned a weight based on the likelihood of the adaptation data. In general, these weights will not be equal, although they can be. To generate a set of weighted runs conditioned on the adaptation data of a participant, follow the instructions described in Performing multiple runs of a COIN model simulation, but in addition pass the adaptation data of the participant to the class object via the 'adaptation' property (the data should be in vector form with one element per trial—use NaN on trials where adaptation was not measured) and set the parameters of the COIN model to their maximum-likelihood (fitted) values via the appropriate properties.
To view view code showing how to fit the parameters of the COIN model to the average data of a group of participants, see fit_parameters_to_data.m (specify the number of participants in the group at the top of the file; BADS installation required). Note that this method assumes that the average adaptation data across participants can be calculated without reordering the adaptation data of each participant (the nth average adaptation data point of the group is calculated as the average nth adaptation data point of each participant).
To fit the parameters of the COIN model to data from multiple experiments simultaneously (e.g., spontaneous and evoked recovery), the negative log-likelihood of each experiment should be calculated and the sum of the log-likelihoods passed to an optimiser.