Models

Model object

class pyddm.model.Model(drift=DriftConstant(drift=0), noise=NoiseConstant(noise=1), bound=BoundConstant(B=1), IC=ICPointSourceCenter(), overlay=OverlayNone(), name='', dx=0.005, dt=0.005, T_dur=2.0, fitresult=None, choice_names=('correct', 'error'))[source]

Bases: object

A full simulation of a single DDM-style model.

Each model simulation depends on five key components:

  • A description of how drift rate (drift) changes throughout the simulation.
  • A description of how variability (noise) changes throughout the simulation.
  • A description of how the boundary changes throughout the simulation.
  • Starting conditions for the model
  • Specific details of a task which cause dynamic changes in the model (e.g. a stimulus intensity change)

This class manages these, and also provides the affiliated services, such as analytical or numerical simulations of the resulting reaction time distribution.

IC(conditions)[source]

The initial distribution at t=0.

Returns a length N ndarray (where N is the size of x_domain()) which should sum to 1.

can_solve_cn(conditions={})[source]

Check whether this model is compatible with Crank-Nicolson solver.

All bound functions which do not depend on time are compatible.

can_solve_explicit(conditions={})[source]

Check explicit method stability criterion

flux(x, t, conditions)[source]

The flux across the boundary at position x at time t.

get_dependence(name)[source]

Return the dependence object given by the string name.

get_fit_result()[source]

Returns a FitResult object describing how the model was fit.

Returns the FitResult object describing the last time this model was fit to data, including the loss function, fitting method, and the loss function value. If the model was never fit to data, this will return FitResultEmpty.

get_model_parameter_names()[source]

Get an ordered list of the names of all parameters in the model.

Returns the name of each model parameter. The ordering is arbitrary, but is uaranteed to be in the same order as get_model_parameters() and set_model_parameters(). If multiple parameters refer to the same “Fittable” object, then that object will only be listed once, however the names of the parameters will be separated by a “/” character.

get_model_parameters()[source]

Get an ordered list of all model parameters.

Returns a list of each model parameter which can be varied during a fitting procedure. The ordering is arbitrary but is guaranteed to be in the same order as get_model_parameter_names() and set_model_parameters(). If multiple parameters refer to the same “Fittable” object, then that object will only be listed once.

get_model_type()[source]

Return a dictionary which fully specifies the class of the five key model components.

has_analytical_solution()[source]

Is it possible to find an analytic solution for this model?

parameters()[source]

Return all parameters in the model

This will return a dictionary of dictionaries. The keys of this dictionary will be “drift”, “noise”, “bound”, “IC”, and “overlay”. The values will be dictionaries, each containing the parameters used by these. Note that this includes both fixed parameters and Fittable parameters. If a parameter is fittable, it will return either a Fittable object or a Fitted object in place of the parameter, depending on whether or not the model has been fit to data yet.

set_model_parameters(params)[source]

Set the parameters of the model from an ordered list.

Takes as an argument a list of parameters in the same order as those from get_model_parameters(). Sets the associated parameters as a “Fitted” object. If multiple parameters refer to the same “Fittable” object, then that object will only be listed once.

simulate_trial(conditions={}, cutoff=True, rk4=True, seed=0)[source]

Simulate the decision variable for one trial.

Given conditions conditions, this function will simulate the decision variable for a single trial. It will cut off the simulation when the decision variable crosses the boundary unless cutoff is set to False. By default, Runge-Kutta is used to simulate the trial, however if rk4 is set to False, the less efficient Euler’s method is used instead. This returns a trajectory of the simulated trial over time as a numpy array.

Note that this will return the same trajectory on each run unless the random seed seed is varied.

Also note that you shouldn’t normally need to use this function. To simulate an entire probability distributions, call Model.solve() and the results of the simulation will be in the returned Solution object. This is only useful for finding individual trajectories instead of the probability distribution as a whole.

simulated_solution(conditions={}, size=1000, rk4=True, seed=0)[source]

Simulate individual trials to obtain a distribution.

Given conditions conditions and the number size of trials to simulate, this will run the function “simulate_trial” size times, and use the result to find a histogram analogous to solve. Returns a Sample object.

Note that in practice you should never need to use this function. This function uses an outdated method to simulate the model and should be used for comparison perposes only. To produce a probability density function of boundary crosses, use Model.solve(). To sample from the probability distribution (e.g. for finding confidence intervals for limited amounts of data), call Model.solve() and then use the Solution.resample() function of the resulting Solution.

solve(conditions={}, return_evolution=False, force_python=False)[source]

Solve the model using an analytic solution if possible, and a numeric solution if not.

First, it tries to use Crank-Nicolson as the solver, and then backward Euler. See documentation of Model.solve_numerical() for more information.

The return_evolution argument should be set to True if you need to use the Solution.get_evolution() function from the returned Solution.

Return a Solution object describing the joint PDF distribution of reaction times.

solve_analytical(conditions={}, force_python=False)[source]

Solve the model with an analytic solution, if possible.

Analytic solutions are only possible in a select number of special cases; in particular, it works for simple DDM and for linearly collapsing bounds and arbitrary single-point initial conditions. (See Anderson (1960) for implementation details.) For most reasonably complex models, the method will fail. Check whether a solution is possible with has_analytic_solution().

If successful, this returns a Solution object describing the joint PDF. If unsuccessful, this will raise an exception.

solve_numerical(method='cn', conditions={}, return_evolution=False, force_python=False)[source]

Solve the DDM model numerically.

Use method to solve the DDM. method can either be “explicit”, “implicit”, or “cn” (for Crank-Nicolson). This is the core DDM solver of this library.

Crank-Nicolson is the default and works for any model with constant bounds.

Implicit is the fallback method. It should work well in most cases and is generally stable.

Normally, the explicit method should not be used. Also note the stability criteria for explicit method is:

noise^2/2 * dt/dx^2 < 1/2

It returns a Solution object describing the joint PDF. This method should not fail for any model type.

return_evolution (default=False) governs whether or not the function returns the full evolution of the pdf as part of the Solution object. This only works with methods “explicit” or “implicit”, not with “cn”.

force_python makes PyDDM use the solver written in Python instead of the optimized solver written in C.

solve_numerical_c(conditions={})[source]

Solve the DDM model using the implicit method with C extensions.

This function should give near identical results to solve_numerical_implicit. However, it uses compiled C code instead of Python code to do so, which should make it much (10-100x) faster.

This does not current work with non-Gaussian diffusion matrices (a currently undocumented feature).

solve_numerical_cn(conditions={})[source]

Solve the DDM model numerically using Crank-Nicolson.

This uses the Crank Nicolson method to solve the DDM at each timepoint. Results are then compiled together. This is the core DDM solver of this library.

It returns a Solution object describing the joint PDF.

solve_numerical_explicit(conditions={}, **kwargs)[source]

Solve the model using the explicit method (Forward Euler).

See documentation for the solve_numerical method.

solve_numerical_implicit(conditions={}, **kwargs)[source]

Solve the model using the implicit method (Backward Euler).

See documentation for the solve_numerical method.

t_domain()[source]

A list of all of the timepoints over which the joint PDF will be defined (increments of dt from 0 to T_dur).

x_domain(conditions, t=None)[source]

A list which spans from the lower boundary to the upper boundary by increments of dx.

Fittable object

class pyddm.model.Fittable(**kwargs)[source]

Bases: float

For parameters that should be adjusted when fitting a model to data.

Each Fittable object does not need any parameters, however several parameters may improve the ability to fit the model. In particular, maxval and minval ensure we do not choose an invalid parameter value. default is the value to start with when fitting; if it is not given, it will be selected at random.

default()[source]

Choose a default value.

This chooses a value for the Fittable object randomly abiding by any constraints. Note that calling this function multiple times will give different results.

Fitted object

class pyddm.model.Fitted(val, **kwargs)[source]

Bases: pyddm.model.Fittable

Parameters which used to be Fittable but now hold a value.

This extends float so that the parameters for the Fittable can be saved in the final model, even though there are no more Fittable objects in the final model.

Sample object (for RT data)

class pyddm.sample.Sample(choice_upper, choice_lower, undecided=0, choice_names=('correct', 'error'), **kwargs)[source]

Bases: object

Describes a sample from some (empirical or simulated) distribution.

Similarly to Solution, this is a glorified container for three items: a list of reaction times for for the two choices (corresponding to upper and lower DDM boundaries), and the number of undecided trials. Each can have different properties associated with it, known as “conditions” elsewhere in this codebase. This is to specifiy the experimental parameters of the trial, to allow fitting of stimuli by (for example) color or intensity.

To specify conditions, pass a keyword argument to the constructor. The name should be the name of the property, and the value should be a tuple of length two or three. The first element of the tuple should be a list of length equal to the number of correct trials, and the second should be equal to the number of error trials. If there are any undecided trials, the third argument should contain a list of length equal to undecided.

By default, the choice associated with the upper boundary is “correct responses” and the lower boundary is “error responses”. To change these, set the choice_names argument to be a tuple containing two strings, with the names of the boundaries. So the default is (“correct”, “error”), but could be anything, e.g. (“left”, “right”), (“high value” and “low value”), etc. This is sometimes referred to as “accuracy coding” and “stimulus coding”. When fitting data, this must match the choice names of the model.

Optionally, additional data can be associated with each independent data point. These should be passed as keyword arguments, where the keyword name is the property and the value is a tuple. The tuple should have either two or three elements: the first two should be lists of properties for the correct and error reaction times, where the properties correspond to reaction times in the correct or error lists. Optionally, a third list of length equal to the number of undecided trials gives a list of conditions for these trials. If multiple properties are passed as keyword arguments, the ordering of the undecided properties (in addition to those of the correct and error distributions) will correspond to one another.

cdf(choice, dt=0.01, T_dur=2)[source]

An estimate of the cumulative density function of sample RTs for a given choice.

choice should be the name of the choice for which to obtain the cdf, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

Note that the return value will not converge to one, but both choices plus the undecided distribution will collectively converge to one.

cdf_corr(dt=0.01, T_dur=2)[source]

The correct component of the joint CDF.

This method is deprecated, use Sample.cdf() instead.

cdf_err(dt=0.01, T_dur=2)[source]

The error (incorrect) component of the joint CDF.

This method is deprecated, use Sample.cdf() instead.

condition_combinations(required_conditions=None)[source]

Get all values for set conditions and return every combination of them.

Since PDFs of solved models in general depend on all of the conditions, this returns a list of dictionaries. The keys of each dictionary are the names of conditions, and the value is a particular value held by at least one element in the sample. Each list contains all possible combinations of condition values.

If required_conditions is iterable, only the conditions with names found within required_conditions will be included.

condition_names()[source]

The names of conditions which hold some non-zero value in this sample.

condition_values(cond)[source]

The values of a condition that have at least one element in the sample.

cond is the name of the condition from which to get the observed values. Returns a list of these values.

static from_numpy_array(data, column_names, choice_names=('correct', 'error'))[source]

Generate a Sample object from a numpy array.

data should be an n x m array (n rows, m columns) where m>=2. The first column should be the response times, and the second column should be the choice that the trial corresponds to. E.g., by default, the choices are (1 == correct, 0 == error), but these can be changed by passing in a tuple of strings to the choice_names variable. E.g. (“left”, “right”) means that (1 == left, 0 == right).

Any remaining columns in data after the first two should be conditions. column_names should be a list of length m of strings indicating the names of the conditions. The order of the names should correspond to the order of the columns. This function does not yet work with undecided trials.

static from_pandas_dataframe(df, rt_column_name, choice_column_name=None, choice_names=('correct', 'error'), correct_column_name=None)[source]

Generate a Sample object from a pandas dataframe.

df should be a pandas array. rt_column_name and choice_column_name should be strings, and df should contain columns by these names.

The column with the name rt_column_name should be the response times, and the column with the name choice_column_name should be the choice that the trial corresponds to. E.g., by default, the choices are (1 == correct, 0 == error), but these can be changed by passing in a tuple of strings to the choice_names variable. E.g. (“left”, “right”) means that (1 == left, 0 == right).

Any remaining columns besides these two should be conditions. This function does not yet work with undecided trials.

correct_column_name is deprecated and included only for backward compatibility.

items(choice=None, correct=None)[source]

Iterate through the reaction times.

choice is whether to iterate through RTs corresponding to the upper or lower boundary, given as the name of the choice, e.g. “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

correct is a deprecated parameter for backward compatibility, please use choice instead.

For each iteration, a two-tuple is returned. The first element is the reaction time, the second is a dictionary containing the conditions associated with that reaction time.

If you just want the list of RTs, you can directly iterate through “sample.corr” and “sample.err”.

mean_decision_time()[source]

The mean decision time in the correct trials.

pdf(choice, dt=0.01, T_dur=2)[source]

An estimate of the probability density function of sample RTs for a given choice.

choice should be the name of the choice for which to obtain the pdf, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

Note that the return value will not sum to one, but both choices plus the undecided distribution will collectively sum to one.

pdf_corr(dt=0.01, T_dur=2)[source]

The correct component of the joint PDF.

This method is deprecated, use Sample.pdf() instead.

pdf_err(dt=0.01, T_dur=2)[source]

The error (incorrect) component of the joint PDF.

This method is deprecated, use Sample.pdf() instead.

prob(choice)[source]

Probability of a given choice response.

choice should be the name of the choice for which to obtain the probability, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s

prob_correct()[source]

The probability of selecting the right response.

This method is deprecated, use Sample.prob() instead.

prob_correct_forced()[source]

The probability of selecting the correct response if a response is forced.

This method is deprecated, use Sample.prob_forced() instead.

prob_error()[source]

The probability of selecting the incorrect (error) response.

This method is deprecated, use Sample.prob() instead.

prob_error_forced()[source]

The probability of selecting the incorrect response if a response is forced.

This method is deprecated, use Sample.prob_forced() instead.

prob_forced(choice)[source]

Probability of a given response if a response is forced.

choice should be the name of the choice for which to obtain the probability, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s

If a trajectory is undecided, then a response is selected randomly.

prob_undecided()[source]

The probability of selecting neither response (undecided).

subset(**kwargs)[source]

Subset the data by filtering based on specified properties.

Each keyword argument should be the name of a property. These keyword arguments may have one of three values:

  • A list: For each element in the returned subset, the specified property is in this list of values.
  • A function: For each element in the returned subset, the specified property causes the function to evaluate to True.
  • Anything else: Each element in the returned subset must have this value for the specified property.

Return a sample object representing the filtered sample.

static t_domain(dt=0.01, T_dur=2)[source]

The times that corresponds with pdf/cdf_corr/err parameters (their support).

to_pandas_dataframe(rt_column_name='RT', choice_column_name='choice', drop_undecided=False, correct_column_name=None)[source]

Convert the sample to a Pandas dataframe.

rt_column_name is the column label for the response time, and choice_column_name is the column label for the choice (corresponding to the upper or lower boundary).

Because undecided trials do not have an RT or choice, they are cannot be added to the data frame. To ignore them, thereby creating a dataframe which is smaller than the sample, set drop_undecided to True.

Solution object (for theoretical RT histograms)

class pyddm.solution.Solution(pdf_choice_upper, pdf_choice_lower, model, conditions, pdf_undec=None, pdf_evolution=None)[source]

Bases: object

Describes the result of an analytic or numerical DDM run.

This is a glorified container for a joint pdf, between the response options (correct, error, and undecided) and the response time distribution associated with each. It stores a copy of the response time distribution for both the correct case and the incorrect case, and the rest of the properties can be calculated from there.

It also stores a full deep copy of the model used to simulate it. This is most important for storing, e.g. the dt and the name associated with the simulation, but it is also good to keep the whole object as a full record describing the simulation, so that the full parametrization of every run is recorded. Note that this may increase memory requirements when many simulations are run.

cdf(choice)[source]

The cumulative density function of model RTs for a given choice.

choice should be the name of the choice for which to obtain the cdf, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

Note that the return value will not converge to one, but both choices plus the undecided distribution will collectively converge to one.

cdf_corr()[source]

The correct component of the joint CDF.

This method is deprecated, use Solution.cdf() instead.

cdf_err()[source]

The error (incorrect) component of the joint CDF.

This method is deprecated, use Solution.cdf() instead.

evaluate(rt, choice=None, correct=None)[source]

Evaluate the pdf at a given response time.

rt is a time, greater than zero, at which to evaluate the pdf.

choice is whether to evaluate on the upper or lower boundary, given as the name of the choice, e.g. “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

correct is a deprecated parameter for backward compatibility, please use choice instead.

Returns the value of the pdf at the given RT. Note that, despite being from a probability distribution, this may be greater than 0, since it is a continuous probability distribution.

mean_decision_time()[source]

The mean decision time in the correct trials (excluding undecided trials).

pdf(choice)[source]

The probability density function of model RTs for a given choice.

choice should be the name of the choice for which to obtain the pdf, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

Note that the return value will not sum to one, but both choices plus the undecided distribution will collectively sum to one.

pdf_corr()[source]

The correct component of the joint PDF.

This method is deprecated, use Solution.pdf() instead.

pdf_err()[source]

The error (incorrect) component of the joint PDF.

This method is deprecated, use Solution.pdf() instead.

pdf_evolution()[source]

The evolving state of the simulation: An array of size x_domain() x t_domain() whose columns contain the cross-sectional pdf for every time step.

If the model contains overlays, this represents the evolving state of the simulation before the overlays are applied. This is because overlays do not specify what to do with the diffusion locations corresponding to undercided probabilities. Additionally, all of the necessary information may not be stored, such as the case with a non-decision time overlay.

This means that in the case of models with a non-decision time t_nd, this gives the evolving probability at time T_dur + t_nd.

If no overlays are in the model, then sum(pdf_corr()[0:t]*dt) + sum(pdf_err()[0:t]*dt) + sum(pdf_evolution()[:,t]*dx) should always equal 1 (plus or minus floating point errors).

Note that this function will fail if the solution was not generated to contain information about the evolution of the pdf. This is not enabled by default, as it causes substantial memory overhead. To enable this, see the documentation for the Model.solve() argument “return_evolution”, which should be set to True.

pdf_undec()[source]

The final state of the simulation, same size as x_domain().

If the model contains overlays, this represents the final state of the simulation before the overlays are applied. This is because overlays do not specify what to do with the diffusion locations corresponding to undercided probabilities. Additionally, all of the necessary information may not be stored, such as the case with a non-decision time overlay.

This means that in the case of models with a non-decision time t_nd, this gives the undecided probability at time T_dur + t_nd.

If no overlays are in the model, then pdf_corr() + pdf_err() + pdf_undec() should always equal 1 (plus or minus floating point errors).

prob(choice)[source]

Probability of a given choice response within the time limit.

choice should be the name of the choice for which to obtain the probability, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s

prob_correct()[source]

Probability of correct response within the time limit.

This method is deprecated, use Solution.prob() instead.

prob_correct_forced()[source]

Probability of correct response if a response is forced.

Forced responses are selected randomly.

This method is deprecated, use Solution.prob_forced() instead.

prob_correct_sign()[source]

Probability of correct response if a response is forced.

Forced responses are selected by the position of the decision variable at the end of the time limit T_dur.

This is only available for the implicit method.

This method is deprecated, use Solution.prob_sign() instead.

prob_error()[source]

Probability of incorrect (error) response within the time limit.

This method is deprecated, use Solution.prob() instead.

prob_error_forced()[source]

Probability of incorrect response if a response is forced.

Forced responses are selected randomly.

This method is deprecated, use Solution.prob_forced() instead.

prob_error_sign()[source]

Probability of incorrect response if a response is forced.

Forced responses are selected by the position of the decision variable at the end of the time limit T_dur.

This is only available for the implicit method.

This method is deprecated, use Solution.prob_sign() instead.

prob_forced(choice)[source]

Probability of a given response if a response is forced.

choice should be the name of the choice for which to obtain the probability, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s

If a trajectory is undecided at the time limit (T_dur), then a response is selected randomly.

prob_sign(choice)[source]

Probability of a given response if a response is forced.

choice should be the name of the choice for which to obtain the probability, corresponding to the upper or lower boundary crossings. E.g., “correct”, “error”, or the choice names specified in the model’s choice_names parameter.

If a trajectory is undecided at the time limit (T_dur), then a response is sampled from the distribution of decision variables at the final timepoint.

prob_undecided()[source]

The probability of not responding during the time limit.

resample(k=1, seed=None)[source]

Generate a list of reaction times sampled from the PDF.

k is the number of TRIALS, not the number of samples. Since we are only showing the distribution from the correct trials, we guarantee that, for an identical seed, the sum of the two return values will be less than k. If no undecided trials exist, the sum of return values will be equal to k.

This relies on the assumption that reaction time cannot be less than 0.

seed specifies the random seed to use in sampling. If unspecified, it does not set a random seed.

Returns a Sample object representing the distribution.