[R] R and MatLab implementations of the same model differs
peter dalgaard
pdalgd at gmail.com
Fri Jul 5 01:38:16 CEST 2013
On Jul 5, 2013, at 00:41 , Jannetta Steyn wrote:
> Hi All
>
> Thanks for all the comments. I have looked at everything you have pointed
> out and this is the situation at the moment:
>
>> 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.
>
> You are right. I have to admit when I referred to equations I referred to
> those in the ST function, and not the initialisation and parameter
> equations. I fixed the line, but look you said it still doesn't produce the
> same results as MatLab.
>
>> 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.
>
> The correct equation is:
>
> iLeak_axon = gLeak_axon_AB.*(v_axon_AB-ELeak_axon_AB);
>
> I think this error was the result of a variable name change at some point
> to make the to scripts the same. I have fixed the MatLab script but it is
> not making any difference.
>
>>
>> Final question for the OP: if the model is supposed to be dynamic, isn't
> it suspicious that Matlab gives a constant result?
>
> Now this is true ... I'll chat to my MatLab expert about this. I was hoping
> that it is not a constant but that it has small variations in the value
> that I can't see on the graph. However, I did not know how to check the
> values ... another question for my MatLab expert tomorrow.
>
>> Also, in this block
>>
>> gNa_axon_AB=300e-3
>> gK_axon_AB=52.5-3
>> gLeak_axon_AB=0.0018e-3
>>
>> the middle line looks like a typo for 52.5e-3 and the 3rd line looks very
> low -- googling suggests values like (120, 36, 0.3) mS/cm^3, which would be
> more consistent with a value around 1e-3.
>
> Yes you are right. 52.5-3 should be 52.5e-3. I made the change but still no
> luck.
>
> With regards to the value of gLeak_axon_AB, the values I'm using for the
> model are from a paper. I tried changing it to what you suggested, but that
> seems to result in a constant value of -60.00006
But isn't -60 = ELeak_axon_AB exactly the right steady-state value for the system with both the K and Na channels closed?
>
> More late night staring ahead for me ...
>
> Regards
> Jannetta
>
>
> --
>
> ===================================
> Web site: http://www.jannetta.com
> Email: jannetta at henning.org
> ===================================
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> 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
More information about the R-help
mailing list