<jupyter_start><jupyter_text>Hodgkin Huxley Model Analysis Description: Hodking Huxley model is a ODE model describing excitation of cells through conductance, voltage and current of sodium and potassium, and input and leaking potentials. The model is a system of 4 ODEs: $$ I = C_m \frac{d}{dt} V_m + g_k n^{4}(V_m - E_K) + g_{Na} m^{3} h (V_m - E_{na}) + g_{l} (V_m - E_l) $$ $$ \frac{dn}{dt} = \alpha_{n}(V_m)(1 - n) - \beta_{n}(V_m)n $$$$ \frac{dm}{dt} = \alpha_{m}(V_m)(1 - m) - \beta_{m}(V_m)m $$$$ \frac{dh}{dt} = \alpha_{h}(V_m)(1 - h) - \beta_{h}(V_m)m $$Equations from: http://www.math.pitt.edu/~bdoiron/assets/ermentrout-and-terman-ch-1.pdf<jupyter_code>## import modules 
import numpy as np
import scipy as sp
from scipy.integrate import odeint
import matplotlib.pyplot as plt<jupyter_output><empty_output><jupyter_text>Setting constants<jupyter_code># Membrane capacitance and potential
C_m = 1.0 
V_m = -70.0

# Max sodium conductance and membrane potential
g_Na = 120.0
V_Na = 56.0

# Max potassium conductance and membrane potential
g_K = 36.0
V_K = -77.0 

# Max leak conductance and membrane potential
g_l = 0.3 
V_l = -60

# Time to integrate over
time = np.arange(0.0, 1000.0, 0.1)

# ion variables 
m = 0.05
h = 0.6
n = 0.32<jupyter_output><empty_output><jupyter_text>Defining model functions<jupyter_code># Input current
def I_inj(time):
    step_current = np.piecewise(time, [(time< 500)*(time>=200), (time< 950)*(time>=700)], [10,50])
    return step_current

# potassium rate functions
def a_n(V_m):
    return 0.01 * (V_m + 55) / (1 - np.exp(-(V_m + 55.0) / 10)) # + 1.E-7)

def B_n(V_m):
    return 0.125 * np.exp(-(V_m + 65.0) / 80)

# sodium rate functions
def a_m(V_m): 
    return 0.1 * (40 + V_m) / (1.0 - np.exp(-(V_m + 40.0) / 10.0)) # + 1.E-7)

def a_h(V_m):
    return 0.07 * np.exp(-(V_m + 65.0) / 20.0)

def B_m(V_m):
    return 4 * np.exp(-(V_m + 65.0) / 18.0)

def B_h(V_m):
    return 1 / (np.exp(-(35 + V_m) / 10) + 1)

def K_current(g_K,n,V_m,V_K):
    return g_K * n**4 * (V_m - V_K)

def Na_current(g_Na, m, h, V_m, V_Na):
    return g_Na * m**3 * h * (V_m - V_Na)
    
def leak_current(g_l, V_m, V_l):
    return g_l * (V_m - V_l)
inj = I_inj(450)
inj<jupyter_output><empty_output><jupyter_text>Setting up model system<jupyter_code>def hh_model(y, t, parameters):
    '''Takes in y dynamics variables, V_m, n, m and h. Array of time steps to integrate over. 
    List of constant parameters, reversal potentials and conductances 

    Returns values for derivatives w.r.t time for Voltage and ion rate coefficience n, m and h. 
    '''
    V_m, n, m, h = y
    g_K, V_K, g_Na, V_Na, g_l, C_m = parameters
    
    # Total current through the membrane    
    dVdt = I_inj(t) - K_current(g_K,n,V_m,V_K) - Na_current(g_Na, m, h, V_m, V_Na) - leak_current(g_l, V_m, V_l) / C_m
 
    # Derivative of n, potassium channel activation, w.r.t. time
    dndt = a_n(V_m) * (1 - n) - B_n(V_m) * n

    # Derivative of m, sodium channel activion, w.r.t. time
    dmdt = a_m(V_m) * (1 - m) - B_m(V_m) * m

    # Derivative of h, sodium channel in-activion, w.r.t. time

    dhdt = a_h(V_m) * (1 - h) - B_h(V_m) * h
    
    return [dVdt, dndt, dmdt, dhdt] 
y0 = [V_m, n, m, h]
parameters = [g_K, V_K, g_Na, V_Na, g_l, C_m]

X = hh_model(y0, time, parameters)
sol = odeint(hh_model,y0, time, args=(parameters,))
Voltage = sol[:,0]
iK = K_current(g_K,n = sol[:,1],V_m = sol[:,0],V_K = V_K)
iNa = Na_current(g_Na, m = sol[:,2], h = sol[:,3], V_m = sol[:,0], V_Na = V_Na)
il = leak_current(g_l, V_m = sol[:,3], V_l = V_l)
plt.plot(time,Voltage)
plt.rcParams["figure.figsize"] = (20,7)
plt.show()
plt.plot(time,sol[:,1])
plt.plot(time,sol[:,2])
plt.plot(time,sol[:,3])
plt.rcParams["figure.figsize"] = (20,7)
plt.show()
plt.plot(time,iK)
plt.plot(time,iNa)
plt.rcParams["figure.figsize"] = (20,7)
plt.show()
plt.plot(time,il)
plt.rcParams['figure.figsize'] = (20,7)
plt.show()
plt.plot(sol[:,1],sol[:,2])


from mpl_toolkits import mplot3d

import numpy as np
import matplotlib.pyplot as plt


fig = plt.figure()
ax = plt.axes(projection="3d")

z_line = sol[:,0]
x_line = sol[:,1]
y_line = sol[:,2]
ax.plot3D(x_line, y_line, z_line, 'blue')

plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.17292161520190025, 'fraction_numerical': 0.03824228028503563, 'mean_word_length': 3.7050279329608937, 'pattern_counts': {'":': 0, '<': 19, '<?xml version=': 0, '>': 19, 'https://': 0, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '84885', 'license': 'MIT', 'n_tokens_mistral': 2028, 'n_tokens_neox': 1891, 'n_words': 456, 'path': 'python_models/Hodgkin_Huxley.ipynb', 'repo_name': 'PavanChaggar/neuronal_models', 'size': '4210'}