Dopamine-modulated medium spiny neuron, reduced model (Humphries et al. 2009)

 Download zip file 
Help downloading and running models
We extended Izhikevich's reduced model of the striatal medium spiny neuron (MSN) to account for dopaminergic modulation of its intrinsic ion channels and synaptic inputs. We tuned our D1 and D2 receptor MSN models using data from a recent (Moyer et al, 2007) large-scale compartmental model. Our new models capture the input-output relationships for both current injection and spiking input with remarkable accuracy, despite the order of magnitude decrease in system size. They also capture the paired pulse facilitation shown by MSNs. Our dopamine models predict that synaptic effects dominate intrinsic effects for all levels of D1 and D2 receptor activation. Our analytical work on these models predicts that the MSN is never bistable. Nonetheless, these MSN models can produce a spontaneously bimodal membrane potential similar to that recently observed in vitro following application of NMDA agonists. We demonstrate that this bimodality is created by modelling the agonist effects as slow, irregular and massive jumps in NMDA conductance and, rather than a form of bistability, is due to the voltage-dependent blockade of NMDA receptors
1 . Humphries MD, Lepora N, Wood R, Gurney K (2009) Capturing dopaminergic modulation and bimodal membrane behaviour of striatal medium spiny neurons in accurate, reduced models. Front Comput Neurosci 3:26 [PubMed]
Model Information (Click on a link to find other models with that property)
Model Type: Neuron or other electrically excitable cell;
Brain Region(s)/Organism:
Cell Type(s): Neostriatum medium spiny direct pathway GABA cell;
Gap Junctions:
Receptor(s): D1; D2; GabaA; AMPA; NMDA;
Transmitter(s): Dopamine;
Simulation Environment: MATLAB;
Model Concept(s): Action Potential Initiation; Parameter Fitting; Simplified Models; Parkinson's; Bifurcation;
Implementer(s): Humphries, Mark D [m.d.humphries at];
Search NeuronDB for information about:  Neostriatum medium spiny direct pathway GABA cell; D1; D2; GabaA; AMPA; NMDA; Dopamine;
%%% script to show properties of neuron models...

clear all

dt = 0.1; % time step

% dopamine level
D1 = 0.8;
D2 = 0.8;

% original MSN parameters
C = 50; vr = -80;  vt = -25; 
k = 1; a = 0.01; b = -20; 
c = -55; d = 150; vpeak = 40; 

load fit_model_results_NEWtuning

% MS neuron parameters in saved file
k = izipars(1); a = izipars(2); b = izipars(3); c = izipars(4); vr = izipars(5); vpeak = izipars(6);

% for control... increase a...
% a = 0.02;

% found MS parameters: X = [C,vt,d]
C = X(1); vt =X(2); d = X(3);

% extra DA model parameters in saved file
KIR = XD1(1);    % KIR modifier 
LCA = XD1(2);    % LCA modifier
vrD1 = vr*(1+D1*KIR);
dD1 = d*(1-D1*LCA);

% D2 - intrinsic
alpha = XD2;
kD2 = k*(1-alpha*D2);

% paired-pulse faciliation (see Mahon paper first for times and
% sizes...)
Ton1 = 50; Toff1 = 250;
ISIs = [200:100:1000];
Ion = 400;
I = 0;
ISIstore = 200;

%% set up sim
deltaT = zeros(numel(ISIs),1); deltaTD1 = zeros(numel(ISIs),1); deltaTD2 = zeros(numel(ISIs),1);
for loop=1:numel(ISIs)
    % set up on and off pulses
    Ton = [Ton1 Toff1+ISIs(loop)];
    Toff = [0 Toff1 Toff1+ISIs(loop)+(Toff1-Ton1)];
    T = Toff(end)+200;
    t = dt:dt:T;
    n = length(t); % number of time points
    v = vr*ones(1,n); u=0*v;
    vD1 = vr*ones(1,n); uD1=0*vD1;
    vD2 = vr*ones(1,n); uD2=0*vD2;
    for i = 1:n-1

        if any(i*dt == Ton)
            I = Ion;
        elseif any(i*dt == Toff)
            I = 0;
        v(i+1) = v(i) + dt*(k*(v(i)-vr)*(v(i)-vt)-u(i) + I)/C;
        u(i+1) = u(i) + dt*a*(b*(v(i)-vr)-u(i));

        if v(i+1)>=vpeak
        %%% D1 type
        vD1(i+1) = vD1(i) + dt*(k*(vD1(i)-vrD1)*(vD1(i)-vt)-uD1(i) + I)/C;
        uD1(i+1) = uD1(i) + dt*a*(b*(vD1(i)-vrD1)-uD1(i));
        % spikes?   
        if vD1(i+1)>=vpeak
            vD1(i)=vpeak; vD1(i+1)=c; 
        %--- D2 type                    
        vD2(i+1) = vD2(i) + dt*(kD2*(vD2(i)-vr)*(vD2(i)-vt)-uD2(i) + I)/C;

        uD2(i+1) = uD2(i) + dt*a*(b*(vD2(i)-vr)-uD2(i));
        % spikes?   
        if vD2(i+1)>=vpeak
            vD2(i)=vpeak; vD2(i+1)=c; 

    % spike times for baseline
    spkts = find(v == vpeak)*dt;
    spk1 = spkts(find(spkts < Toff(2),1,'first'));
    spk2 = spkts(find(spkts > Ton(2) & spkts <Toff(end),1,'first'));
    t1 = spk1 - Ton(1); t2 = spk2 - Ton(2);
    if ~isempty(t1)& ~isempty(t2)  deltaT(loop) = t1-t2; end
    % for D1 model
    spktsD1 = find(vD1 == vpeak)*dt;
    spk1 = spktsD1(find(spktsD1 < Toff(2),1,'first'));
    spk2 = spktsD1(find(spktsD1 > Ton(2) & spktsD1 <Toff(end),1,'first'));
    t1 = spk1 - Ton(1); t2 = spk2 - Ton(2);
    if ~isempty(t1) & ~isempty(t2) deltaTD1(loop) = t1-t2; end
    % for D2 model
    spktsD2 = find(vD2 == vpeak)*dt;
    spk1 = spktsD2(find(spktsD2 < Toff(2),1,'first'));
    spk2 = spktsD2(find(spktsD2 > Ton(2) & spktsD2 <Toff(end),1,'first'));
    t1 = spk1 - Ton(1); t2 = spk2 - Ton(2);
    if ~isempty(t1)& ~isempty(t2) deltaTD2(loop) = t1-t2; end
    % draw current pulse on figure
    poff = -95;
    pon = poff + 2*((Ion/1000)/0.1);
    pulse = [ones(Ton(1)/dt,1).*poff; ones((Toff(2)-Ton(1))/dt,1).*pon; ones((Ton(2)-Toff(2))/dt,1).*poff;...
                ones((Toff(3)-Ton(2))/dt,1).*pon; ones((T-Toff(3))/dt,1).*poff];

    figure(loop); clf
    subplot(131), plot(t,v,'k','LineWidth',2) 
    hold on, plot(t,pulse,'LineWidth',2); 
    xlabel('Time (ms)'); ylabel('Membrane potential (mV)')
    subplot(132), plot(t,vD1,'k','LineWidth',2) 
    hold on, plot(t,pulse,'LineWidth',2); 
    xlabel('Time (ms)'); ylabel('Membrane potential (mV)')
    subplot(133), plot(t,vD2,'k','LineWidth',2) 
    hold on, plot(t,pulse,'LineWidth',2); 
    xlabel('Time (ms)'); ylabel('Membrane potential (mV)')
    if ISIs(loop) == ISIstore
        baseoutput = [t' v' pulse];

hold on
% plot D1 and D2 effects

save paired_pulse_responses

Loading data, please wait...