Contexte:
J'utilise le package Ryacas , en essayant de trouver une expression symbolique pour la variance à grand échantillon d'un estimateur MLE .
Pour ce faire, j'ai besoin du (inverse de la) matrice Hessienne de la fonction de vraisemblance. Je n'ai pas de mathématique (et la version en ligne me semble trop lourde à utiliser) et j'essaye donc avec le package Ryacas qui est une interface avec l'ordinateur YACAS système d'algèbre.
Question:
Cependant, je n'arrive pas à comprendre comment calculer le Hessian . En utilisant le guide ici: https: //cran.r -project.org/web/packages/Ryacas/vignettes/high-level.html me donne une erreur pour cela. Voici un exemple minimal reproductible (tiré de ce guide)
L <- yac_symbol("x^2 * (y/4) - a*(3*x + 3*y/2 - 45)")
y_fn(L, "HessianMatrix")
Quand j'essaye d'une autre manière, à savoir en utilisant la (nouvelle?) Interface
L <- yac_symbol("x^2 * (y/4) - a*(3*x + 3*y/2 - 45)")
Hessian(L)
Error in Hessian(L) : could not find function "Hessian"
Je n'arrive pas non plus à obtenir une réponse utilisable, mais juste HessianMatrix ((x ^ 2 * y) / 4-a * (3 * x + (3 * y) / 2-45 ))
Quelqu'un aurait-il une idée de comment résoudre ce problème? Serait grandement apprécié!
Merci
3 Réponses :
L'erreur, Error in Hessian (L): n'a pas pu trouver la fonction "Hessian" , est due au fait que la fonction Hessian ne fait pas partie de Ryacas , il fait à la place partie du package numDeriv dans R.
Vous trouverez également ci-dessous la documentation de la fonction Hessian dans R:
https://www.rdocumentation.org/packages/numDeriv /versions/2016.8-1.1/topics/hessian
J'espère que cela résout votre problème. Veuillez donc d'abord installer le package numDeriv , puis utiliser la fonction Hessian .
J'ai également essayé d'utiliser l'exemple ci-dessus, voici les résultats: p>
{}
Sortie:
library(Ryacas)
library(numDeriv)
L <- yac_symbol("x^2 * (y/4) - a*(3*x + 3*y/2 - 45)")
Hessian(L)
La sortie est un ensemble vide, mais je suppose que vous pouvez vous référer à la documentation de Hessian et cela aidera certainement.
Merci, mais j'ai bien peur que ce ne soit pas ce que je recherche. Voici pourquoi: si vous faites défiler vers le bas du lien vers la documentation Ryacas, vous verrez la fonction "Hessian". De plus, le lien vers la fonction que vous donnez montre que la fonction est nommée différemment (avec un petit h et pas un grand), et, enfin, la fonction hessian du package numDeriv donne une approximation numérique - pas ce que je recherche ici ( comme le titre l'indique, je recherche un calcul symbolique).
library(Ryacas)
yac_str("HessianMatrix(x^2 * (y/4) - a*(3*x + 3*y/2 - 45), {x,y})")
# "{{y/2,x/2},{x/2,0}}"
yac_str("PrettyForm(HessianMatrix(x^2 * (y/4) - a*(3*x + 3*y/2 - 45), {x,y}))") %>% cat
# / \
# | / y \ / x \ |
# | | - | | - | |
# | \ 2 / \ 2 / |
# | |
# | / x \ ( 0 ) |
# | | - | |
# | \ 2 / |
# \ /
Ref: yacas.readthedocs
Merci, Stéphane, cela a fonctionné pour l'exemple du jouet! Malheureusement, lorsque j'utilise ma fonction réelle (c'est un peu plus compliqué), j'obtiens beaucoup d'expressions "Deriv" en Hesse qui ne sont pas évaluées. PrettyForm n'est pas utile dans ce cas car je souhaite traiter davantage la sortie et en soustraire une autre expression pour la variance de grand échantillon d'un autre estimateur (en essayant de montrer qu'un estimateur est plus efficace qu'un autre pour toutes les valeurs du vrai (multinomial) distribution sous-jacente). Merci, je marque votre réponse comme correcte, cela m'a donné quelques conseils sur la façon de procéder!
@Jean_N: Merci, bon à savoir. Je mentionnerai aussi simplement le package caracas pour vous; ceci est basé sur SymPy plutôt que sur yacas.
Voici ce que vous pouvez faire avec Ryacas v1.1.1:
> library(Ryacas)
> packageVersion("Ryacas")
[1] ‘1.1.1’
> x <- ysym("x")
> y <- ysym("y")
> a <- ysym("a")
> L <- x^2 * (y/4) - a*(3*x + 3*y/2 - 45)
> H <- Hessian(L, c("x", "y", "a"))
> H
{{ y/2, x/2, -3},
{ x/2, 0, (-3)/2},
{ -3, (-3)/2, 0}}
> as_r(H)
expression(rbind(c(y/2, x/2, -3), c(x/2, 0, -3/2), c(-3, -3/2,
0)))
> eval(as_r(H), list(x = 2, y = 2, a = 2))
[,1] [,2] [,3]
[1,] 1 1.0 -3.0
[2,] 1 0.0 -1.5
[3,] -3 -1.5 0.0
Je suis revenu sur cette question beaucoup plus tard en recommençant à travailler sur le problème. En attendant, j'ai dû passer à nlsDeriv pour m'aider à résoudre mon problème, mais votre solution fonctionne et est exactement ce dont j'avais besoin, merci! J'espère travailler davantage avec Ryacas maintenant et mieux le connaître.