Page 1 of 1

help: how do you: check value of parameter or record a value

Posted: Fri Jun 15, 2007 11:41 am
by shrbert724
I've been playing around with record() feature but can't quit get it all to work. I want to know the value of certain parameters at a specific time so I can use that value as an initial condition to refine my model.

below is the code that I'm working with. how would i get neuron to output or save the instantaneous value of a variable at a specific point. I'm interested in the values of h, n, r, sout, ina, ik to name a few if that helps. I know this would be pretty simple in other coding languages so I assume it is here too, but when I try and code it, it doesn't work.

//hoc code
create stn
objref app_stn, app_gpe, ge2sn1, sn2ge1, apstn, apgpe, vstn, vgpe
stn {
insert STN

el_STN = -60.0
ek = -80.0
ena = 55.0
eca = 140.0
ett = 140.0
eahp = -80.0


cm = 100
diam = 5.5
L = 11

app_stn=new ISine(0.5)
app_stn.Ibase = 4
app_stn.Iamp = 0
app_stn.Ifreq = 0
app_stn.Ibase2 = 4
app_stn.Iamp2 = 0
app_stn.Ifreq2 = 0
app_stn.Ibase3 = 4
app_stn.Iamp3 = 0
app_stn.Ifreq3 = 0
app_stn.del = 200
app_stn.dur = 700
vstn = new Vector()
apstn = new APCount(0.5)
apstn.record(vstn)
}


//stn mod file
NEURON {
SUFFIX STN
USEION k READ ek WRITE ik
USEION na READ ena WRITE ina
NONSPECIFIC_CURRENT il
USEION ca READ eca WRITE ica
USEION tt READ ett WRITE itt VALENCE 2
USEION ahp READ eahp WRITE iahp VALENCE 1
RANGE gk_bar, gna_bar, gl, gtt_bar, gca_bar, gahp_bar, ik, ina, il, itt, ica, iahp, n_inf, m_inf, h_inf, a_inf, b_inf, r_inf, s_inf
}

UNITS {
(S) = (siemens)
(mV) = (millivolts)
(mA) = (milliamps)
(M) = (moles/liter)
(mM) = (0.001 M)
(mC) = (millicoulombs)
(uC) = (microcoulombs)
}

PARAMETER {
gl = 0.225 (S/cm2)
gk_bar = 4.5 (S/cm2)
gna_bar = 3.75 (S/cm2)
gtt_bar = 0.05 (S/cm2)
gca_bar = 0.05 (S/cm2)
gahp_bar = 0.9 (S/cm2)
ek = -80.0 (mV)
ena = 55.0 (mV)
eca = 140.0 (mV)
el = -60.0 (mV)
ett = 140.0 (mV)
eahp = -80.0 (mV)
theta_m = -30.0 (mV)
theta_h = -39.0 (mV)
theta_n = -32.0 (mV)
theta_r = -67.0 (mV)
theta_a = -63.0 (mV)
:theta_b needs to be unitless for y_inf definition
theta_b = 0.4
theta_s = -39.0 (mV)
theta_g = 30.0 (mV)
sigma_m = 15.0 (mV)
sigma_h = -3.1 (mV)
sigma_n = 8.0 (mV)
sigma_r = -2.0 (mV)
sigma_a = 7.8 (mV)
:sigma_b needs to be unitless for y_inf function
sigma_b = -0.1
sigma_s = 8.0 (mV)
k1 = 15.0 (mM)
kca = 22.5 (C/M s cm2)
epsilon = 3.75e-5 (M cm2/mC)
: connection terms (A_alpha, B_alpha in Terman p.232)
a_alpha = 5.0 (/ms)
b_alpha = 1.0 (/ms)
:from paper, theta_heavy changed to 10 because could be either of two values from paper, so a rough average taken as 10 mV
theta_heavy = 10.0 (mV)
sigma_heavy = 8.0 (mV)

: %%% Theta_Tau %%%%%%%%%%%
theta_tau_h = -57.0 (mV)
theta_tau_n = -80.0 (mV)
theta_tau_r = 68.0 (mV)

: %%% Sigma_Tau %%%%%%%%%%%%
sigma_tau_h = -3.0 (mV)
sigma_tau_n = -26.0 (mV)
sigma_tau_r = -2.2 (mV)

: %%% Tau1 (msec) %%%%%%%%%%%
tau1_h = 500.0 (ms)
tau1_n = 100.0 (ms)
tau1_r = 17.5 (ms)

: %%% Tau0 (msec) %%%%%%%%%%
tau0_h = 1.0 (ms)
tau0_n = 1.0 (ms)
tau0_r = 40.0 (ms)

: % Phi (unitless) %%%%%%%%%%%%%%%
phi_h = 0.75
phi_n = 0.75
phi_r = 0.2
}


ASSIGNED {
v (mV)
ik (mA/cm2)
ina (mA/cm2)
ica (mA/cm2)
itt (mA/cm2)
il (mA/cm2)
iahp (mA/cm2)
n_inf m_inf h_inf a_inf b_inf r_inf s_inf
}




STATE { h n r sout siout cai (mM) }

BREAKPOINT {
SOLVE states METHOD cnexp
ss()
ik = gk_bar * n^4 * (v - ek)
il = gl * (v - el)
ina = gna_bar * m_inf^3 * h * (v - ena)
itt = gtt_bar * a_inf^3 * b_inf^2 * (v - ett)
ica = gca_bar * s_inf^2 * (v - eca)
iahp = gahp_bar * (v - eahp) * (cai / (cai + k1))
}



INITIAL {
ss()
h = h_inf
n = n_inf
r = r_inf
cai = 5.5
sout = 0
}



DERIVATIVE states {
n' = phi_n * (n_inf - n) / taun(tau0_n, tau1_n, v, theta_tau_n, sigma_tau_n)
h' = phi_h * (h_inf - h) / tauh(tau0_h, tau1_h, v, theta_tau_h, sigma_tau_h)
r' = phi_r * (r_inf - r) / taur(tau0_r, tau1_r, v, theta_tau_r, sigma_tau_r)
cai' = epsilon * (-ica - itt - kca * cai)
sout' = a_alpha * (1 - sout) * heavyside(v, theta_heavy, sigma_heavy) - b_alpha * sout
}



PROCEDURE ss() {
m_inf = x_inf(v,theta_m,sigma_m)
h_inf = x_inf(v,theta_h,sigma_h)
n_inf = x_inf(v,theta_n,sigma_n)
r_inf = x_inf(v,theta_r,sigma_r)
a_inf = x_inf(v,theta_a,sigma_a)
s_inf = x_inf(v,theta_s,sigma_s)
b_inf = y_inf(r,theta_b,sigma_b)
}

FUNCTION heavyside(Vm (mV), theta_heavy (mV), sigma_heavy (mV)) {
heavyside = 1/(1+exp(-(Vm-theta_heavy)/sigma_heavy))
}

FUNCTION tauh(tau0_h (ms), tau1_h (ms), Vm (mV), theta_tau_h (mV), sigma_tau_h (mV)) (ms) {
:HELP sigma_tau_h/n are both negative while there is a negative in front of the Vm statement
tauh = tau0_h+(tau1_h/(1+exp(-(Vm-theta_tau_h)/(sigma_tau_h))))
}

FUNCTION taun(tau0_n (ms), tau1_h (ms), Vm (mV), theta_tau_n (mV), sigma_tau_n (mV)) (ms) {
taun = tau0_n+(tau1_n/(1+exp(-(Vm-theta_tau_n)/(sigma_tau_n))))
}

FUNCTION taur(tau0_r (ms), tau1_r (ms), Vm (mV), theta_tau_r (mV), sigma_tau_r (mV)) (ms) {
taur = tau0_r+(tau1_r/(1+exp(-(Vm-theta_tau_r)/(sigma_tau_r))))
}

FUNCTION x_inf(Vm (mV), theta (mV), sigma (mV)) {
x_inf = 1 / (1 + exp(-(Vm - theta)/(sigma)))
}

FUNCTION y_inf(r, theta_b, sigma_b) {
y_inf = 1 / (1 + exp((r - theta_b)/(sigma_b))) - 1 / (1 + exp(-(theta_b/sigma_b)))
}

Re: help: how do you: check value of parameter or record a v

Posted: Fri Jun 15, 2007 8:57 pm
by ted
shrbert724 wrote:I've been playing around with record() feature but can't quit get it all to work.
Too much code. The general algorithm for doing something new is:
1. Read the documentation first.
2. Try it on a simple example until you get it to work.
3. Only after it works in a simple setting should one try it in a more complex situation.

If you want to record the time course of a variable:
The Programmer's Reference entry on the Vector class's record() method is fairly explicit,
and includes an example--
http://www.neuron.yale.edu/neuron/stati ... tml#record
Try it with a single compartment model that contains only the hh or even just the pas mechanism. Use minimal code until you are sure it works.