R and MatLab implementations of the same model differs
On 04-07-2013, at 19:56, peter dalgaard <pdalgd at gmail.com> wrote:
On Jul 4, 2013, at 19:11 , Berend Hasselman wrote:
On 04-07-2013, at 18:42, Jannetta Steyn <jannetta at henning.org> wrote:
Hi Ben and others I don't quite know how to explain the "doesn't work" in more detail without any visual aid.
You said that R got into an indefinite loop, whatever that maybe.
When you run the two scripts it is easy to see the difference. MatLab produces a line on x= -55. This is what I expect - a more or less straight line. R on the other hand the result drops from -55 (the initial value) to -80 and then goes up to -71.37092. The two scripts have exactly the same equations.
I don't think so. In the R script you have init = c(v_axon_AB=-55,mNa_axon_AB=1,hNa_axon_AB=0,mK_axon_AB=1) That is not the same as in your Matlab script. To make them the same you should replace the line with init = c(v_axon_AB=-55,mNa_axon_AB=0,hNa_axon_AB=1,mK_axon_AB=0) Using this line gives quite different results. But not the same as Matlab. It seems that you ought to carefully check all your equations. I don't know enough about Matlab/Octave syntax to determine if the results of the two systems are identical.
Also, the line iLeak_axon_AB <- gLeak_axon_AB*(v_axon_AB-ELeak_axon_AB) differs from the Matlab code. Changing it gave me all-NA results, though. And the with() construct assumes that init and parms are named vectors; parms is not, so the values are taken from the global environment. It is not clear whether that makes a difference.
Same NA results here. Making parms a named vector makes no difference but is is certainly necessary to do that. In that case there is probably an inconsistency in the file simulate.m near the end: The last non empty line before the line out = ? reads iLeak_axon = gLeak_axon_AB.*(v_axon_AB-ELeak_axon_AB); and this doesn't agree with the line iLeak_axon = ELeak_axon_AB*(v_axon_AB-ELeak_axon_AB); in the function xdot. Final question for the OP: if the model is supposed to be dynamic, isn't it suspicious that Matlab gives a constant result? Berend
Berend
I have even named the variables the same in the two scripts and copied the equations across to make sure I haven't made any typos. Can one attach images to posts? I'll try. The flatline image is the plot from MatLab and the other is the plot from R. Thanks Jannetta On 4 July 2013 16:52, Ben Bolker <bbolker at gmail.com> wrote:
Berend Hasselman <bhh <at> xs4all.nl> writes:
On 04-07-2013, at 17:15, Jannetta Steyn <jannetta <at> henning.org>
wrote:
Hi folks I have implemented a model of a neuron using Hodgkin Huxley equations
in
both R and MatLab. My preference is to work with R but R is not giving
me
the correct results. I also can't use ode45 as it just seems to go
into an
indefinite loop. However, the MatLab implementation work fine with the same equations, parameters and initial values and ode45. Below is my R and MatLab implementations.
No problem in running your R file. Have plot. (Mac mini Core2Duo 2.66Ghz running R 3.0.1 Patched (2013-06-19 r62992) on Mac OS X 10.8.4: 16.5 seconds.) Trying to run your Matlab file in Octave. No success yet because of unavailable ode45.
I'm impressed that you (BH) went to the trouble of checking on this vague "doesn't work" question. If you want to go farther I think you can get ode45 for octave by installing the odepkg package: http://octave.sourceforge.net/odepkg/index.html
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
-- =================================== Web site: http://www.jannetta.com Email: jannetta at henning.org =================================== <Rplot01.png><MatPlot01.png>______________________________________________ R-help at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
______________________________________________ R-help at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
-- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Email: pd.mes at cbs.dk Priv: PDalgd at gmail.com