Neurons are computers: The Computational Power of Dendrites

By: Tom Breloff

Re-posted from: http://www.breloff.com/Neurons-are-computers/

The human brain is composed of billions of computational units, called neurons. Neurons transmit information to each other through electrical pulses, or spikes. A common misrepresentation is that the neuron is the basic conputational building block; that spikes are transferred and integrated at the soma, and that axons and dendrites are simply a network of electrical cables connecting the neurons together. Treating the dendrite as cabling is convenient; the mathematics are simpler, and it’s easier to reason about the role of neurons, both individually and in groups. As we’ll see, this simplification is wrong, and it glosses over a critical computational component in our neural architecture.

If you missed it, please check out my first post: Efficiency is Key: Lessons from the Human Brain

Some of the graphics and animations below were created in Julia with Plots.jl. Check out the IJulia notebook.

Figure A: Components of a Neuron [1]

The Basics of Spiking

A neuron’s output is a sequence of electrical bursts; termed action potentials, or spikes. A neuron’s cell body (soma) accumulates post-synaptic (after the synaptic cleft) chemical energy in the form of a membrane potential. At the site of the axon hilcock, net input potentials are summated, and a spike is generated. This spike travels quickly through short sections of the axon surrounded by myelin sheath, where the signal is then boosted at each node (the gaps between sheaths). This successive boosting allows the current to quickly reach the tips of the axonal branches.

Axons (output) form connections with the dendrites (input) of other neurons; this connection is called a synapse. When a signal reaches the end of an axon at a synapse, it triggers the release of neurotransmitters which then cross the synaptic gap, binding and opening channels on the dendrite. The transfer of chemicals (for example sodium and potassium) are what change the interal charge of the dendrite, and also the soma. This charge builds up until crossing a theoretical membrane threshold and triggering a run-away process which results in a spike. (See this short video for a visual explanation of the process)

Axons terminate at many pre-synaptic (before the synaptic cleft) sites. Synaptic terminals are possible in many configurations. Axons can synapse directly with the soma (axosomatic), leading to direct potential transfer, or they could even synapse with another synaptic terminal (axosynaptic), leading to modulatory effects. However the most common connection is axon-to-dendrite (axodendritic).


Figure B: Synaptic Connection Locations [1]

Location, location, location…

The location of the synaptic terminal has a large effect on the resulting change in membrane potential (which, when high enough, will produce a spike). Those synapses connecting axon to soma (axosomatic) have nearly direct energy transfer; the postsynaptic response (i.e. how the receiving neuron responds) will be somewhat linear when compared with the presynaptic input.

For synapses on the dendritic tree, the postsynaptic potential (PSP) is affected both by passive properties of dendrites [2], and synaptic inhibition from neighboring synapses [3], both driving a sublinear response. This means that postsynapitic potentials have diminishing marginal response as more synapses on the dendritic segment are active.

In addition, dendritic synapses have active properties. When synapses are active in isolation, the relative response at the soma is generally small. However, when multiple synapses become active with spatio-temporal proximity (meaning nearby in both space and time) then it can trigger a supralinear response: a dendritic spike. This active conductance is relatively large as compared to the individial postsynaptic responses; however it quickly saturates, and additional presynaptic input has little effect.


Figure C: Membrane Potential Response Shapes

This means that, depending on the structure and properties of a neuron’s dendritic tree and synaptic connectivity, membrane potential can be effected in highly nonlinear ways. Notice in Figure C how the sublinear line is the square root of input, while the supralinear line is similar to the sigmoid function. This is intentional to emphasize parallels to current deep learning methods. For more details on the sublinear and supralinear curves, see [3].

Dendritic Segments as a Key Computational Component

As hinted before, nearby synapses interact in nonlinear ways when they are active at the same time. It turns out that relative distance is not the defining characteristic; synapses are grouped into dendritic segments which are somewhat isolated from each other. Depending on the properties of the segment, as well as the location, the synaptic interaction can be linear, sublinear, or supralinear.

Remember also that dendrites have a branching, tree-like structure. One advantage is that the branch points, where the tree splits, have inhibitory properties. Thus the dendritic forks can integrate or inhibit potentials from more distal dentritic areas.

To summarize: A neuron has dendritic branches protruding from the soma. Those branches split at dendritic forks into more branches. On each branch are some number (possibly zero) of dendritic segments. A segment may contain many axodendritic synapses. Integration, inhibition, and amplification can occur at every level of the heirarchy, creating a rich and deep computation structure within a single neuron. In Julia:

abstract NeuralComponent

type Synapse <: NeuralComponent
    presynaptic_neuron
    presynaptic_current
end

type Segment <: NeuralComponent
    synapses
end

type Branch <: NeuralComponent
    segments
    fork
end    

type Soma
    inputs::Set{NeuralComponent}
    membrane_potential::Float64
end

To be continued… (when I have time)

[1] Blausen.com staff. “Blausen gallery 2014”. Wikiversity Journal of Medicine. doi: 10.15347/wjm/2014.010. ISSN 20018762.

[2] Abrahamsson, T, Cathala, L, Matsui, K, Shigemoto, R, and Digregorio, D. A. (2012). Thin dendrites of cerebellar interneurons confer sublinear synaptic integration and a gradient of short-term plasticity. Neuron. 73, 1159-1172. doi: 10.1016/j.neuron.2012.01.027

[3] Tran-Van-Minh A, Caze RD, Abrahamsson T, Cathala L, Gutkin BS and DiGregorio DA (2015) Contribution of sublinear and supralinear dendritic integration to neuronal computations. Front. Cell. Neurosci. 9:67. doi: 10.3389/fncel.2015.00067

Using Julia’s C Interface to Utilize C Libraries

By: Christopher Rackauckas

Re-posted from: http://www.stochasticlifestyle.com/using-julias-c-interface-utilize-c-libraries/

I recently ran into the problem that Julias’s (GNU Scientific Library) GSL.jl library is too new, i.e. some portions don’t work as of early 2016. However, to solve my problem I needed access to adaptive Monte Carlo integration methods. This means it was time to go in depth in Julia’s C interface. I will go step by step on how I created and compiled the C code, and called it from Julia.

What I wished to use was the GSL Vegas functions. Luckily, there is a good example on this page for how to use the library. I started with this code. I then modified it a bit. To start, I changed the main function header to:

int monte_carlo_integrate (double (*integrand)(double *,size_t,void *),double* retRes,double* retErr,int mode,int dim,double* xl,double* xu,size_t calls){

Before it took in no values, but I will want to be able to do most things from Julia. First of all, I changed the function name from main to avoid namespace issues. This is just good practice when working with shared libraries. then I added a bunch of arguments. The first argument is a function pointer, where the function is defined just as g is in the example. Then I pass in pointers which we become the return values. I add mode (1,2,3 are the different integration types), dim, xl, xu, and calls as passed in values to easily extend the functionality. The rest of the C-code stays mostly the same: I change the parts which reference g to integrand (I like the function name better) and comment out the print functions (they tend to no print correctly). So I end up with a C-function as follows:

#include <stdlib.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_monte.h>
#include <gsl/gsl_monte_plain.h>
#include <gsl/gsl_monte_miser.h>
#include <gsl/gsl_monte_vegas.h>
int monte_carlo_integrate (double (*integrand)(double *,size_t,void *),double* retRes,double* retErr,int mode,int dim,double* xl,double* xu,size_t calls){
  double res, err;
  const gsl_rng_type *T;
  gsl_rng *r;
  gsl_monte_function G = { *integrand, dim, 0 };
  gsl_rng_env_setup ();
  T = gsl_rng_default;
  r = gsl_rng_alloc (T);
  if (mode==1){
    gsl_monte_plain_state *s = gsl_monte_plain_alloc (dim);
    gsl_monte_plain_integrate (&G, xl, xu, dim, calls, r, s,
                               &res, &err);
    gsl_monte_plain_free (s);
    /* display_results ("plain", res, err); */
  }
  if (mode==2){
    gsl_monte_miser_state *s = gsl_monte_miser_alloc (dim);
    gsl_monte_miser_integrate (&G, xl, xu, dim, calls, r, s,
                               &res, &err);
    gsl_monte_miser_free (s);
    /* display_results ("miser", res, err); */
  }
  if (mode==3){
    gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (dim);
    gsl_monte_vegas_integrate (&G, xl, xu, dim, 10000, r, s,
                               &res, &err);
    /* display_results ("vegas warm-up", res, err); */
    /* printf ("converging...n"); */
    do
      {
        gsl_monte_vegas_integrate (&G, xl, xu, dim, calls/5, r, s,
                                   &res, &err);
        /* printf ("result = % .6f sigma = % .6f "
                "chisq/dof = %.1fn", res, err, gsl_monte_vegas_chisq (s)); */
      }
    while (fabs (gsl_monte_vegas_chisq (s) - 1.0) > 0.5);
    /* display_results ("vegas final", res, err); */
    gsl_monte_vegas_free (s);
  }
  gsl_rng_free (r);
  *retRes = res;
  *retErr = err;
  return 0;
}

Notice how nothing was specifically changed to be “Julia-style”, this is all normal C-code. So how do I use it? First I need to compile it to a shared library. This means I use the -shared tag and save it to a .so. I will need the GSL libraries loaded when compiling, so I add the -lgsl, -lgslcblas, and -lm (math) tags to link those libraries. Then I add the -fPIC tag since Julia’s documentation says so, and tell it the file to compile. The complete code is:

gcc -Wall -shared -o libMonte.so -lgsl -lgslcblas -lm -fPIC monteCarloExample.c

This will give you a .so file. Now we need to use it. I will describe passing in the function last. Here’s the setup. The first of the other variables are the two pointers retRes and retErr where the results will be saved. In Julia, to get a pointer, I make two 1-element arrays as follows:

x = Array{Float64}(1)
y = Array{Float64}(1)

The only other peculiar thing I had to do was to change pi from Julia’s irrational type to a Float64 (Cdouble) and use that to make the array. So for the other variables I used:

fpi = convert(Float64,pi)
xl = [0.; 0.; 0.]
xu = [fpi; fpi; fpi]
calls = 500000
mode = 3
dim = 3

Now I pass this all to C as follows:

ccall((:monte_carlo_integrate,"/path/to/library/libMonte.so"),Int32,(Ptr{Void},Ptr{Cdouble},Ptr{Cdouble},Int32,Int32,Ptr{Cdouble},Ptr{Cdouble},Csize_t),integrand_c,x,y,mode,dim,xl,xu,calls)

The first argument a tuple where the first place is the symbol which says which function in the library to use and the second place is the library name or the path to the library. The second argument is the return type (here it’s Int32 because our function returns int 0 when complete). The next portion is a tuple which defined the types of the arguments we are passing. The first is Ptr{Void} which is used for things like function pointers. Next there are two Cdouble pointers for retRes and retErr. Then two integers, two more pointers, and a Csize_t. Lastly we put in all of the variables we want to pass in.

Recall that the result was stored into the pointers for the second and third variable passed in. These are the pointers x and y. So to get the values of res and err, I de-reference them:

res=x[1]
err=y[1]

Now res and err hold the values for the solution and the error.

I am not done yet since I didn’t talk about the function! To do this, we define the function in Julia. We have to use parametric types or Julia will yet at us, and we have to ensure that the returned value is a Cdouble (in our case). So we re-write the g function in Julia as:

function integrand{T,T2,T3}(x::T,dim::T2,params::T3)
  A = 1.0 / (pi * pi * pi)
  return A / (1.0 - cos(unsafe_load(x,1))*cos(unsafe_load(x,2))*cos(unsafe_load(x,3)))::Cdouble
end

Notice that instead of x[1], we have to use the Julia function unsafe_load(x,1) to de-reference a pointer. However, since this part is in Julia, other things are much safer, like how we can use pi directly without having to convert it to a float. Also notice that we can add print statements in this function, and they will print directly to Julia. You can use this to modify the display_results function to be a Julia function which prints. However, this itself is still not able to be passed into C. To do that, you have to translate it to a C function:

integrand_c = cfunction(integrand,Cdouble,(Ptr{Cdouble},Ptr{Cdouble},Ptr{Cdouble}))

Here we used the Julia cfunction function where the first argument is the function, the second argument is what the function returns, and the third argument is a tuple of C-types that the function will take on. If you look back at the ccall, this integrand_c is what we passed as the first argument to the actual C-function.

Check to see that this all works. It worked for me. For completeness I will put the full Julia code below. Happy programming!

x = Array{Float64}(1)
y = Array{Float64}(1)
fpi = convert(Float64,pi)
xl = [0.; 0.; 0.]
xu = [fpi; fpi; fpi]
 
function integrand{T,T2,T3}(x::T,dim::T2,params::T3)
  A = 1.0 / (pi * pi * pi)
  return 3*A / (1.0 - cos(unsafe_load(x,1))*cos(unsafe_load(x,2))*cos(unsafe_load(x,3)))::Cdouble
end
 
integrand_c = cfunction(integrand,Cdouble,(Ptr{Cdouble},Ptr{Cdouble},Ptr{Cdouble}))
 
calls = 500000
mode = 3
dim = 3
ccall((:monte_carlo_integrate,"/home/crackauc/Public/libMonte.so"),Int32,(Ptr{Void},Ptr{Cdouble},Ptr{Cdouble},Int32,Int32,Ptr{Cdouble},Ptr{Cdouble},Csize_t),integrand_c,x,y,mode,dim,xl,xu,calls)
res=x[1]
err=y[1]

The post Using Julia’s C Interface to Utilize C Libraries appeared first on Stochastic Lifestyle.