Coding up Literature (research papers)


I am curious if anyone has any experience converting equations from research papers into functioning models using python. I am currently trying to implement the math described in this paper, but don’t really know where to start.
If anyone has any tips or insight into how to think about this task, please let me know!

Thanks in advance :slight_smile:

There’s a paywall on that, I think you’re more likely to get help with specific questions.

Do you have a screenshot of the math?

Hello both! Thank you for following up on my post.
Here is a link that you should be able to access the PDF through.
And here are the 2 main equations I would like to implement in Python!

Let me know if you have done anything similar or have any suggestions!!

1 Like

That is way over my head, though we could begin with pencil to paper. Help us with the symbols and describe them in some fundamental way that lets us compute with them.

Being so heavily dependent on symbols, the above math might be easier to work with in Pysym (symPy?), a symbolic computation library. Again, this is beyond my level of comprehension, but given generalized concepts, we can handle transposing of that to individual arithmetic steps. We’re obviously not going to solve the above in one line of code. We need to break it down into basic components and according to how they relate, devise the code to solve each bit.

Seems to me a while back I wrote an algo for solving quadratics. Will see if I can track it down. Might be a while, though.


As @mtf mentioned you might have an easier time with SymPy. I think the simplest starting point would be to create a function to mimic the mathematical function(s) provided and then to simply create a range of points equivalent to those on the axes shown and calculate the result of your f(x,y,z,...). Then you can sanity check your result to those provided and you’ll know if you’re on the right track for any further results.

There are a number of constants you can pre-calculate (the paper provides many if not all) and the plotting seems to be based on eq.2 with additional variable values provided in the figure caption. Once these are all replaced in the function it should be a fairly simple equation.

So find out exactly which equation is being plotted. Simplify the function as much as possible with pre-calculated values (you should have a limited number of actual variable values remaining). Then try and re-create the given plots, you can consider the later results and 3D plots once you’ve successfully re-created some of the given functions.


Thank you @mtf and @tgrtim for your helpful suggestions! Also, much of this is over my head too, but I am slowly catching on, lol.

In terms of creating a function to mimic the one provided, in what sense do you mean ‘mimic’. Are you suggesting I define the function (2) with the same variables but then define variables with values based on some points given on the graph?
So something along these lines:

V_tilde = #given graph value

E_in = -2

#channel open state probability 
P_V = #define function based on equation given, but create additional variables above associated with this equation

#main equation function (2)
def main_equation_2(V_tilde, E_in, P_V):
        I = (V_tilde  -  E_in)*(P_V)
        return I

#####further plotting functions for visualization#####

How might I make use of SymPy with what I have written above, given what I have is in the right direction?

1 Like

You’re welcome, @taylorrayne, to whatever morsel we might leave.

Forgive me drilling right down to the basics of a function.

c = a + b

The left side is a function of the right side.

f(a, b) = a + b

c = f(6, 7)

Counting on you to bring this connect to your deconstruction. Keep in mind iteration. You are plotting a sample set.

1 Like

Yes, you can write functions in Python that act like their mathematical counterparts, there’s a lot of overlap. Comp sci goes into in more detail but things like the following should make it quite clear- Pure function - Wikipedia. I think you’ll find this task easier if you do it that way.

For a loose hint on getting started if you follow the author’s description you’ll find that for the first couple of graphs at least there is only one variable. Everything else is provided so re-creating those first plots is quite straightforward.

Obviously the later plots are more complex but the process would probably be similar. For multi-variates you’d probably want to start working with numpy to make everything easier though once again it is not essential.

If you’ve not looked into it before then SymPy provides some quick examples of how functions can be made a little cleaner when using it. It’s not essential but it’s worth a look-in:

Alternatively since it’s Python3 and the default source encoding is UTF-8 you could swap out the lengthy variable names for the matching characters Ĩ, and so on. Normally that’s a bit sketchy but for this particular task I think it might make the code more readable.

1 Like

Once again, thank you both :slight_smile: (@mtf, @tgrtim)

So far, I have gone ahead and drafted the first equation for inward-rectifying channels.
Your comments on iteration and the single variable have been very helpful!
Although I am facing an error when I try to run through my loop and subtract an integer value from the range of V_tilde. Here is my code below. Note that I have not yet tidied anything - so it is very messy… also, I still have to do the plotting. But let me know what you think so far!

import math #defining variable values # membrane potential of inward-rectifying voltage gated channel ## NEED TO CENTRE MID-POINT AT ZERO V_tilde = range(9) # equilibrium potential E_in = -2 # gating effective charge z_in = 3 z_out = -3 # threshold potentials ## NEED TO UNDERSTAND WHY TWO v_tilde_th_1 = 0.5 v_tilde_th_2 = 2 # channel openstate probibilty for inward-rectifying voltage gated channels def P_in_open_equation(V_tilde): P_in_open1 = 1/(1 + math.exp(z_in*(V_tilde - v_tilde_th_1))) P_in_open2 = 1/(1 + math.exp(z_in*(V_tilde - v_tilde_th_2))) return P_in_open1, P_in_open2 # current potential equation # for threshold 1 for v in V_tilde: I_tilde1 = (V_tilde - E_in)*P_in_open1 # for threshold 2 for v in V_tilde: I_tilde2 = (V_tilde - E_in)*P_in_open2

Also, @tgrtim, how do you think the concept of Pure functions applied to the process above?

Seems like you’ve got a decent plan now. I’m not really sure about the benefit of returning multiple values from that function but the implementation is up to you.

As for the errors, I think the complexity of the names may be tripping you up. Double check what you pass and when.

Regarding the functions what do you think? If you used your function to obtain a value for Ĩ does it match up to the result the mathematical function produces? If you called the same function is the output identical? What about condition 1? Does it even matter? It’s a suggestion for what I’ve found easier to at least start with in the past but at the end of the day do whatever works.

Trying to perfect code at the expense of creating working code is something I’ve tried to shake off myself. It’s often much easier to get your program to work as intended and then to clean up the code when you’re doing something new and unfamiliar.


Alright, I am back with an update.

I have been able to replicate the graphs a) and b) from section II. Here they are:

They are far from aesthetically pleasing and my code is not very pretty either, but I agree with what you mention above, @tgrtim, concerning perfection vs functionality.

On the topic of Pure Functions, I understand the premise behind them, but I am just having a hard time imagining any function that would not be pure…

Thank you both for your guidance!! This has been fun for me and has really helped me learn :slight_smile:
(if you are interested, here is a link to my code. If you take a peak, feel free to comment on anything or critique)


Thanks for sharing that article, it was interesting how electronics mimic biology.
I saw you were trying to work out how to do this:
V_tilde = range(9) - 4

If you’re familiar with functions, you might build this as a function for different inputs.

You could also ditch the variable V_tilde and use list comprehensions
It does work as is though. :smiley:

1 Like

It won’t be long in finding one. They abound. Functions have a way of being impure. The simplest way I can find to describe a pure function is the calculus of linear regression. Every term is a pure function. Think on that, and when you have the time we might go down this avenue further.

def f(x): return x

x = 42
y = f(x)
y == x   #  True

The simplest pure function, identity. It scales up from there.

One neglected to add that maths and programming are each a juxtaposition of the other. They share paradigms and logic, but not a lot else. Ergo, impure functions that maths would not entertain.

In the raw, lambda functions are pretty rare, but as returns from factory functions they are common. This is another avenue to go down. Those functions would also be pure since they do not interact with the outside. They are given a parameter and deal with it. Same outcome every time for the same argument.

The factory function could roll into a scenario where an iterator comes forward. It could manufacture the callback function. Another avenue to go down.

Pure functions and iterators give us a way to write expression based code. The language interprets expressions as values. The more we can learn to write expressive code, the easier it will be to give the computer the values we need it to see.

This is technically a pure function…

My solution - #3 by mtf

1 Like

Thanks for looking over the code! Also, I am glad you found interest in the article :wink:

And yes, I was attempting to work with a range instead, it seems more straightforward, rather than creating a list - that just seems bulky, in my opinion.

Only, I couldn’t work out the TypeError I was getting: Unsupported operand type(s) for -: ‘range’ and ‘int’
Do you have any ideas of what might be going on here or how one could think about solving this issue?
Do you think that your idea with list comprehensions could save the day?

Something like:

V_tilde = [ for v in range(9) ...]

Not really sure how to set this up.
Also, I am not sure I follow what exactly ditching the V_tilde variable could look like.
Any suggestions are appreciated!

1 Like

Ah, this is lovely. Thank you for this thorough expansion on the topic. I am just about to begin liner regression as part of a stats course I am taking in python, currently. You have given me a lot of food for thought - much to chew on, I’ll bookmark this post!

1 Like

I guess I got kinda hasty in posting a response. :upside_down_face:
For the list range I was thinking x = [x - 4 for x in range(9)]
It’s a bit more direct to just do this though:
P_in_open1 = [1/(1 + math.exp(z_in*((v-4) - v_tilde_th_1))) for v in range(9)]

You could make a function of it to allow for different temperatures, or to have single point returns:

def I_tilde(eq_pot, gate_charge, v_thresh, v_tilde):
  I_tilde = v_tilde #your function work here
  return I_tilde

t1_in = [I_tilde(-2, 3, 0.5, x - 4) for x in range(9)]
1 Like

Ah, I love how you set this up as a function and with the list comprehension. Very slick. I will try and give this a go when I have the time to look over the code soon here again.


Pay very close attention.


Hi all. If you are still interested, I have updated the code here: Google Colaboratory

Currently working out how to have multiple initial conditions for the ODE that is graphed in figure 6 of the paper - if you have any ideas, let me know!