-
Notifications
You must be signed in to change notification settings - Fork 6
Home
- The COIN class
- Running the model
- Integrating out observation noise
- Storing variables
- Plotting internal representations
- Inferring internal representations of the COIN model fit to adaptation data
- Fitting the model to data
- Parallel computing
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). These default values 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 simulate learning on a spontaneous recovery paradigm, define a series of perturbations (use NaN to indicate a channel trial):
obj.perturbations = [zeros(1,50) ones(1,125) -ones(1,15) NaN(1,150)];
Then run the model by calling the simulate_COIN method on object obj:
S = obj.simulate_COIN;
The adaptation and state feedback (perturbation plus random observation noise) are contained in S.runs{1}. They can be plotted as follows:
figure
plot(S.runs{1}.y,'.')
hold on
plot(S.runs{1}.yHat)
legend('state feedback','adaptation')
The actual observation noise that a participant perceives is unknown. Therefore, rather than performing a single run of a simulation conditioned on one particular noise sequence, multiple runs of the simulation can be performed, each conditioned on a different noise sequence. The observation noise can then be integrated out by averaging inferences across runs. Use the 'runs' property to specify the number of runs to perform. For example, to perform 2 runs:
obj.runs = 2;
S = obj.simulate_COIN;
S.runs is a cell array with one cell per run and S.weights is vector specifying the relative weight of each run. In this example, each run is assigned an equal weight. However, if adaptation data is passed to the object before running the model, each run will be assigned a weight based on how well it explains the data (see Inferring internal representations fit to adaptation data).
Online inference does not require all of the past values of all inferred variables to be stored in memory. Therefore, to reduce memory requirements, the past values of variables are only stored in memory if they are needed for later analysis (the adaptation and state feedback variables are the only exceptions to this, as they are always stored). To store the values of particular variables on all trials, add the names of these variables to the store property. For example, to store the Kalman gains and responsibilities, define store as
obj.store = {'k','cFilt'};
The store property must be set before calling the simulate_COIN method. For a full list of the names of variables that can be stored, call help COIN
(see 'VARIABLES' section).
The stored variables are contained in the cell array S and can be analysed after the model has run. For example, to compute the Kalman gain of the context c* with the highest responsibility:
for run = 1:obj.runs % loop over runs
for trial = 1:numel(obj.perturbations) % loop over trials
for particle = 1:obj.particles % loop over particles
[~,j] = max(S.runs{run}.cFilt(:,particle,trial));
Kalman_gain(trial,run,particle) = S.runs{run}.k(j,particle,trial);
end
end
end
Within each run, a simple average across particles can be computed (as all particles within the same run have the same weight), whereas across runs, a weighted average should be computed:
figure
plot(sum(mean(Kalman_gain,3).*S.weights,2))
ylabel('Kalman gain | c^*')
To plot specific variables or distributions of variables, activate the relevant flags in the properties of the class object. If plotting a continuous distribution, also specify the points at which to evaluate the distribution. For example, to plot the distribution of the state of each context and the predicted probabilities:
obj.plot_state_given_context = true;
obj.state_values = linspace(-1.5,1.5,500); % values to evaluate state | context at
obj.plot_predicted_probabilities = true;
These properties should be set before calling the simulate_COIN method so that the variables needed to generate the plots can be stored. For a full list of variables or distributions that can be plotted and how to plot them, see the 'plot flags' and 'plot inputs' sections of properties.
After the model has run, call the plot_COIN method:
P = obj.plot_COIN(S);
This will generate the requested plots—a state | context plot and a predicted probabilities plot in this example. The structure P contains the data that is plotted (to see how the data in P 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 variables or distributions are averaged across particles and runs. In general, using more runs will result in cleaner, less variable results.
The COIN model can be fit either to the data of an individual participant or to the average data of a group of participants. Here, the group case is considered, as an individual participant is a special case of a group with 1 member.
The COIN model is fit to data by finding the parameters that minimise the negative log of the likelihood function. To calculate the negative log-likelihood, create an object array with one object per participant:
obj(1,8) = COIN; % array of 8 objects for a group of 8 participants
For each object, define the model parameters, the paradigm (perturbations and, if applicable, cues) and the adaptation data of the participant using the properties of the object. The adaptation property should be a vector with one element per trial (use NaN on trials where adaptation was not measured). Note that the number of adaptation measurements should be the same for all participants (the nth average adaptation measurement is calculated as the average nth adaptation measurement across participants).
After the properties of the object array have been set, call the objective_COIN method:
negative_log_likelihood = obj.objective_COIN;
This returns a stochastic estimate of the negative log-likelihood, which can be passed to an optimiser. The estimate is stochastic because it is calculated from simulations that are conditioned on random observation noise. To aid parameter optimisation, the variance of the estimate can be reduced by increasing the number of runs using the 'runs' property (this is best done in conjunction with parallel computing to avoid excessive runtimes). An optimiser that can handle a stochastic objective function should be used (e.g., BADS).
Note that to fit the model to multiple paradigms simultaneously (e.g., spontaneous and evoked recovery), the negative log-likelihood should be calculated for each paradigm separately and their sum passed to the optimiser.
To see how to fit the COIN model to synthetic force-field adaptation data with sensory cues, see fit_model_to_data.m (note that this example requires installation of BADS).
When multiple runs of a simulation are executed using parameters fit to data, each run can be assigned a weight based on how well it explains the data. In general, these weights will not be equal (although they can be, as weights are reset when runs are resampled during particle filtering). To generate a set of weighted runs, set the model parameters to their maximum likelihood estimates, pass the data to the class object via the adaptation property and call the simulate_COIN method. The resultant weights should be used to compute a weighted average of variables across runs (as in Analysing stored variables)
The time it takes to execute multiple runs in series (e.g. using a for loop) can be prohibitively long if the number of runs is large, as computation scales linearly with the number of runs. Therefore, in practice, it is advisable to use the smallest number of runs that achieves desired performance. In addition, runtime can be reduced by executing multiple runs in parallel (e.g., across different CPU cores on a computer cluster). To execute multiple runs in parallel, specify the maximum number of CPU cores available using the maxCores property. The default setting of maxCores is 0, which executes multiple runs in series.
Parallel computing is utilised by the simulate_COIN and objective_COIN methods but not with the plot_COIN method.