[R] eigenvalues of matrices of partial derivatives with ryacas
Adam Zeilinger
zeil0006 at umn.edu
Fri Mar 16 19:10:13 CET 2012
Dear Gabor,
Thank you so much for your help! Your suggested code worked!
I have a couple of questions. I'm trying to understand your code so I
can use it in future programming needs. First, my understanding of the
first line of your code is that it creates an R function "EigenValues"
that corresponds to the Yacas function of the same name:
EigenValues <- function(x) Sym("EigenValues(", x, ")")
The "x" variable is kind of a place holder for the product of the
matrices for which I want the eigenvalues, correct? Can you explain the
use of the two sets of quotes in the Sym function? Why does the x have
to be outside quotes and why did you include the commas on either side
of the x?
Second, the term "xx^2" appears in the result of the eigenvalue call.
If I ignore this term, this answer is equivalent to the answer I get in
Mathematica. How should I interpret this term, since x only appears in
the "EigenValue" R function and is not part of my model?
Finally, I am using TinnR 2.3.6.3 as an R-GUI. When I'm using Ryacas
functions, I can send single lines of code to R via TinnR, but I get the
following error message when I send multiple lines of code:
> source(.trPaths[4], echo=TRUE, max.deparse.length=150)
Error in srcfilecopy(filename, lines, file.info(filename)[1, "mtime"]) :
unused argument(s) (file.info(filename)[1, "mtime"])
Any thoughts on what the problem might be and how to troubleshoot it?
Thanks again for your help!
Adam
On 3/15/2012 5:16 PM, Gabor Grothendieck wrote:
> On Thu, Mar 15, 2012 at 7:51 PM, Adam Zeilinger<zeil0006 at umn.edu> wrote:
>> Hello,
>>
>> I am trying to construct two matrices, F and V, composed of partial
>> derivatives and then find the eigenvalues of F*Inverse(V). I have the
>> following equations in ryacas notation:
>>
>>> library(Ryacas)
>>> FIh<- Expr("betah*Sh*Iv")
>>> FIv<- Expr("betav*Sv*Ih")
>>> VIh<- Expr("(muh + gamma)*Ih")
>>> VIv<- Expr("muv*Iv")
>> I successfully found the partial derivatives:
>>
>>> f11<- deriv(FIh, "Ih")
>>> f12<- deriv(FIh, "Iv")
>>> f21<- deriv(FIv, "Ih")
>>> f22<- deriv(FIv, "Iv")
>>> v11<- deriv(VIh, "Ih")
>>> v12<- deriv(VIh, "Iv")
>>> v21<- deriv(VIv, "Ih")
>>> v22<- deriv(VIv, "Iv")
>> Next I would like to put these partial derivatives into two matrices, F and
>> V:
>>
>>> F<- Expr("{{f11, f12}, {f21, f22}}")
>>> V<- Expr("{{v11, v12}, {v21, v22}}")
>> Finally, I would like to find the eigenvalues of F*Inverse(V). Something
>> like:
>>
>>> yacas("EigenValues(F*Inverse(V))")
>> However, this does not work. I get the following error message:
>>
>> In function "While" : bad argument number 1 (counting from 1)The offending
>> argument $ii49<= $nr49 evaluated to Not Length-1<0CommandLine(1) : Invalid
>> argument
>>
>> According to Mathematica, the correct eigenvalues are:
>>
>> {-((Sqrt[betah] Sqrt[betav] Sqrt[Sh] Sqrt[Sv])/Sqrt[gamma muv + muh muv]),
>> (Sqrt[betah] Sqrt[betav] Sqrt[Sh] Sqrt[Sv])/Sqrt[gamma muv + muh muv]}
>>
>> I don't understand the error message. Any suggestions on how to get the
>> correct eigenvalues using R would be greatly appreciated. I'm using R
>> 2.14.0.
> Try doing it this way instead:
>
> library(Ryacas)
> EigenValues<- function(x) Sym("EigenValues(", x, ")")
>
> Iv<- Sym("Iv"); Ih<- Sym("Ih")
> Sv<- Sym("Sv"); Sh<- Sym("Sh")
> betav<- Sym("betav"); betah<- Sym("betah")
> muv<- Sym("muv"); muh<- Sym("muh")
> gamma<- Sym("gamma")
>
> FIh<- betah*Sh*Iv
> FIv<- betav*Sv*Ih
> VIh<- (muh + gamma)*Ih
> VIv<- muv*Iv
>
> f11<- deriv(FIh, Ih)
> f12<- deriv(FIh, Iv)
> f21<- deriv(FIv, Ih)
> f22<- deriv(FIv, Iv)
> v11<- deriv(VIh, Ih)
> v12<- deriv(VIh, Iv)
> v21<- deriv(VIv, Ih)
> v22<- deriv(VIv, Iv)
>
> F<- List(List(f11, f12), List(f21, f22))
> V<- List(List(v11, v12), List(v21, v22))
>
> EigenValues(F * Inverse(V))
>
> For the last line I get:
>
>> EigenValues(F*Inverse(V))
> expression(Roots(xx^2 - betav * Sv * muv * (betah * Sh * (muh +
> gamma))/((muh + gamma) * muv)^2))
>
--
Adam Zeilinger
Post Doctoral Scholar
Department of Entomology
University of California Riverside
www.linkedin.com/in/adamzeilinger
More information about the R-help
mailing list