Numerical solver design

Numerical solver design is a topic that I know next to nothing about. I always say I know a tiny bit of MATLAB very well, but lately in my push to become MATLAB-independent, I’ve been trying to remove bits of proprietary code from my own simulations. (My transition from MATLAB to Octave failed in part because of memory issues in Octave. So I’m now looking at Python as my potential free software savior.)

The biggest hurdle has been my usage of the suite of numerical solvers for ordinary differential equations (ODEs) in MATLAB. These are fast, vectorizable to some extent, and take advantage of some parallel processing on multicore machines. However, using their solvers presents a number of problems when one wants to understand the ins and outs of their algorithms. It’s not that this information isn’t there (using ‘type’ or ‘help’), but the problem with going through code that isn’t yours is an ugly one that isn’t going away soon. Additionally, it’s never been completely clear under what circumstances MATLAB’s solvers start to make use of multiple cores and why that process should be accurate, despite seeing evidence of multicore usage at some level of scale in my own code, with ode23s. And there are some things that I don’t know if MATLAB’s solvers can do well. For instance, I have MxNxT matrices of a dynamical variable connecting M model cells of one type to N model cells of another type, saved over all T time steps. Tracking all of these variables and the some ~20 other variables of varying lengths becomes a real code nightmare in vector forms without name labels, which are available in structs. Of course, there is speed in vectors that isn’t available in structs — speed that I am willing to sacrifice for sanity.

The real bottom line for me is being able to write some MATLAB code that I can easily translate into a Python environment eventually. This process will take a lot longer than I had hoped. By writing solvers now, I hope to be able to take them with me and get some insight into solver design for future projects.

Using structs in a solver. In MATLAB’s ODE solvers, they require a strict type of vectorization that has time in one dimension and all of the dynamical variables in the other dimension. I have synaptic weight matrices that are MxNxT, for M of cell A going to N of cell B, throughout time T. Because of the dimensions of these are tricky to code in MATLAB directly without a lot of transposes and other confusing situations, I was hoping to use structs to name my matrices and make it easy to write vector field equations in a clear manner without having to write confusing converters. It’s not completely out of laziness that I want to avoid doing this; rather, another conversion step introduces one more place where error can be introduced.

It turns out that using structs is not a good idea with the separate solver/vector field design, since it seems like you’d spend too much computational time translating anyway. Translating once in a difficult way but gaining speed is far better than translating thousands of times in an “easy” way, at some unknown expense of speed. I might try next a two-matrix approach: one for class A cells and one for class B cells, which does not make this a general solver but will allow for easier vectorization.

Numerical ODE Solvers

In MATLAB, the suite of numerical solvers for systems of ordinary differential equations are functions that accept variables, including other functions. Since I learned to program MATLAB for my numerical solutions, I am accustomed to this system of passing a function of equations and variables to the ODE solver function. While there is something gratifying in my mind about the distinct separation and, in theory, plug-and-play ease of this system, I cannot help but think that it’s computationally more efficient to have the ODE solver integrated into the numerical application. When people say that they write their own ODE solvers, I always assumed that they meant an ODE solver that accepted a function file of equations and iterated over them, in much the same way as MATLAB. However, it’s finally becoming clear to me that this is not the case, and by ODE solvers, most mean that they are simply hard-coding a numerical algorithm or scheme for iterating their equations. The difference is subtle and may boil down to a difference in computational versus mathematical approaches. In my current project, I am interested in interacting rhythms in biophysical models of neural systems, which requires me to solve large systems of ODEs numerically.

While I make the transition from MATLAB to GNU Octave/R, I’m learning about what implementations are reasonable to translate and what new things I need to learn. Of course, I am interested in maximizing efficiency in the shift, but I have to be careful not to do so at the expense of lost computational efficiency, which happened in my first attempt at using ‘lsode’ as a drop-in solver replacement for a complicated system of ODEs. Since I have decided to ditch graphics in MATLAB in lieu of another system that has fewer four letter words associated with it, the ODE solver was my last true obstacle in the translation from MATLAB to Octave. Now that I am going to solve that problem by “writing my own” scheme, I am that much closer to being completely free and open source in my science.

Evolved to run?

In a new book called Born to Run, author Christopher McDougall explores an interesting theory that humans have evolved to be distance runners. Woven into the story of a running tribe called the raramuri of the Copper Canyons of Mexico, McDougall explores the research from human evolutionary biologists Daniel Lieberman and Dennis Bramble. The basic theory is that humans have evolved as bipedal animals with the specific capacity for endurance running that is not shared by other primates such as chimpanzees. It’s a very seductive idea to distance runners to think that we’re just doing what we so naturally do. In fact, we have evolved to do exactly this, so we’re really doing it old school! And to top it off, the McDougall book has gotten hundreds of runners interested in barefoot running, which is put forth as a far more natural way of running than in overbuilt shoes. Biomechanically this is a seductive idea as well, one that I’ve had fair success in implementing in my own running. (I’m up to around 12 miles comfortably in nothing but foot-form fitting slippers.)

Recently the find of an Ardipithecus ramidus, or Ardi, was unveiled to the public. Ardi is 4.4 million years old (Myo), far older than Lucy (an Australopithecus afarensis), who is a mere 3.3 Myo. I wonder whether or not the features of Homo sapiens (us!) that are posited to be beneficial for endurance running are found in these precursors to Homo.

For science and running geeks interested in reading some of the work done on the evolution of humans as endurance runners, see this Nature review by Lieberman and Bramble. A few other articles by the same authors exist, accessible via PubMed, and that review points to several good places to look for other relevant work.

Psychosomatic health

I can’t qualify posts like this enough; they are zeroth order thoughts that are more like armchair musings than they are even close to anything scientific.

Fairly often I come across a science news article about some study that correlates behavioral/mental health with some kind of somatic manifestation. There’s a term that many bandy around, “psychosomatic,” often used superciliously to indicate that one’s pains are somehow made up in one’s mind. However, pain is a very personal, neurological experience that makes objectivity difficult to assess. A recent headline I saw posits a link between loneliness in women and a higher incidence of breast cancer. To assume without qualification, for a moment, that there is a mechanism behind this observation, I was thinking about what that might be. The intuitive reaction is often to assume that it may be through inaction that the brain/body might not fight aggressive cancer cell growth if the loneliness is associated with less than optimal brain function (as is easy to imagine). However, what happens when we think about this in the opposite way? What if the mechanism of action is similar to programmed organism death? What if, in detecting a sub-optimal neurological state, the brain actively contributes to a condition that makes it susceptible to a parasitic toxin such as cancer cells?

This is all firmly in the domain of a gedanken experiment at the moment, but I bet there is some research that investigates several of these issues and potentially how they might work together to understand this. To extend these vagaries further, understanding may lead to effective treatment that could be as simple as “being social” or “making oneself happy.” Mechanism unknown, perhaps, these kinds of simple things may have far greater ramifications for our psychosomatic health (in the non-pejorative way).

Update: the original article will be published in a journal called Cancer Prevention Research.

Meta comments – and some rude ones

It’s difficult to find time to write with some reasonable estimation of quality on a variety of topics in which I’m interested. I’ve meant to write about the Tour de France, the confirmation process of Sonia Sotomayor, and a variety of scientific topics. One difficulty I’ve found is not being able to discuss my work in a meaningful manner on an informal website such as this. I prefer to let my formal talks and (hopefully soon) published papers to communicate interesting results or explored ideas, and I think that peer reviewed science is still the best model, despite something so seductively democratic as the arXiv model of publication, which would be a nightmare to sort through. So there’s an apparent lack of science-related stuff on here, simply for these reasons.

One thing I think I can write about, briefly, is the publication process that I’m currently learning about for the first time in my scientific career. I recently submitted a manuscript with two collaborators to a journal that specializes in computational science, and it was out for peer review about a month ago. We recently received comments back from just the first of two external reviewers, who are not affiliated with the journal directly but presumably experts in the area of research that we are in.

For this first reviewer, the respectfulness of the written comments and questions was greatly appreciated, and I proceeded to incorporate changes and answer questions from that reviewer. The second reviewer apparently did not submit comments to the authors (us), but upon inquiring about this, it turns out that the second reviewer simply did not fill out the review form correctly and the comments to the editors were allowed to be forwarded. The tone of the comments was more caustic, almost accusatory, though there were several valuable points that were made. The reviewer made a few comments that were not at all constructive, simply providing an opinion without any possible suggestions for meeting his or her otherwise arbitrary criteria. I refrain from giving examples here due to confidentiality. Finally, the reviewer was simply wrong about the guidelines of the journal itself in at least one comment that was made, which makes satisfactory revision on this point impossible. Given the tone of this reviewer’s argument, I expect a difficult task in appeasing him or her for the revised submission.

Cocktail party concerto

Zeroth order thought.

Someone recently was listening to a Brahms’ Violin Concerto and asked the question of whether or not the soloist could be picked out amid the background of the much louder orchestra, when playing tutti. I recently had seen some evidence presented by a speaker that might apply here, and I thought it would be fun to consider the two things together.

Essentially, the problem is the following. You have a very large contingent of instruments playing, with a single instrument playing something different. We would like to hear the soloist, either in harmony with the other parts or distinctly, above the din. Now, there is clearly an absolute volume at which the soloist would be drowned out — we believe this from experience. However, assuming we’re not at that level yet, there’s sort of an interesting problem of differentiating a very loud section of instruments that is often of equal timbre but of higher amplitude.

The cocktail party problem refers to the observation that, when in a crowded room in which people are speaking, by looking at a single speaker, one can often hear what she is saying, despite the distracting, surrounding noise. So the concerto problem is similar in this regard.

Dream job post

I was looking at job offerings and came across a post for a neurobiophysicist. Because it contained both the roots “neuro” and “biophysics,” my interest was piqued, since I’ve never heard of anyone referred to by this compounded title. It turns out that I’m apparently training to be a neurobiophysicist, and despite the location of the position (San Antonio, Texas), I am very intrigued by the job post. My hope is that other academic places are looking for someone who fits this description (emphasis added):

We seek to recruit a vigorous academic with an established track-record of research in the broad area of the physics of neurobiological processes. This might be anyone from a “patch clamper” to a computational biophysicist to someone studying signal transduction mechanisms but not necessarily limited to these examples. Trinity University has secured funds to endow the position and provide continuing funding for the position after our HHMI award expires. This professorship will be at the level of full professor. The position comes with a dedicated half-time support staff position, discretionary funds, and a reduced teaching load (6 rather than 9 contact hours per semester).

This position will be an appointment in the Department of Physics, providing some support to introductory physics courses, a course in biophysics (suitable for both physics and neuroscience majors), and a course in the faculty member’s specific area of expertise (specifically supporting neuroscience majors) in the teaching repertoire. In addition, the person hired in the position will maintain a vigorous research program that will engage undergraduate student researchers. The neuroscience major requires independent research as a culminating experience for the degree. Further details about the position are available in the position announcement.

The position will also work with faculty from biology, chemistry, and psychology in supporting the neuroscience major. This interdisciplinary major was established in 2005 with funding from our 2004 HHMI award. Through a combination and retirement and this HHMI-funded neurobiophysicist position, the neuroscience program will be expanded and redefined. In addition to this position, we will be hiring a cognitive scientist in Psychology and an animal behaviorist in Biology to support the neuroscience program.

It’s so accurately describes what I am interested in (and doing) because of aspects of combining computational approaches with electrophysiological techniques. Additionally, the appointment is in a physics department (my undergraduate degree is in physics), with teaching duties in both neuroscience and physics. I’m pretty far from doing a job search, but here’s to hoping positions like this become more common by the time I graduate!