It's pretty hard to debug something when all we are told is "it is
giving me an error".
Take a specific example, (e.g. which we know will be dead flat at (0,0),
just as a check). Here's one way to do it:
> fxy <- quote(x^2*log(x^2+y^2+1) + sin(x+y))
> fxy
x^2 * log(x^2 + y^2 + 1) + sin(x + y)
> Hxy <- function(x, y) NULL
> body(Hxy) <- call('cbind',
fxx = D(D(fxy, "x"), "x"),
fxy = D(D(fxy, "x"), "y"),
fyy = D(D(fxy, "y"), "y"))
> Hxy
function (x, y)
cbind(fxx = 2 * log(x^2 + y^2 + 1) + 2 * x * (2 * x/(x^2 + y^2 +
1)) + (2 * x * (2 * x/(x^2 + y^2 + 1)) + x^2 * (2/(x^2 +
y^2 + 1) - 2 * x * (2 * x)/(x^2 + y^2 + 1)^2)) - sin(x +
y), fxy = 2 * x * (2 * y/(x^2 + y^2 + 1)) - x^2 * (2 * x *
(2 * y)/(x^2 + y^2 + 1)^2) - sin(x + y), fyy = x^2 * (2/(x^2 +
y^2 + 1) - 2 * y * (2 * y)/(x^2 + y^2 + 1)^2) - sin(x + y))
> xy <- expand.grid(x = (-5):5, y = (-2):2)
> cbind(xy, do.call(Hxy, xy))
x y fxx fxy fyy
1 -5 -2 13.0149369 0.87920882 1.87920882
2 -4 -2 11.1066815 0.08339629 0.66389516
3 -3 -2 9.0947006 -0.34667938 -0.40790387
4 -2 -2 7.2919676 0.23085183 -0.65803706
5 -1 -2 5.2801945 1.25223112 0.03000890
6 0 -2 4.1281733 0.90929743 0.90929743
7 1 -2 5.9805455 -0.26964013 0.73035987
8 2 -2 8.0487701 -0.98765432 0.09876543
9 3 -2 9.2121539 -1.45371588 -0.29045058
10 4 -2 10.4767996 -1.27210922 0.03401323
11 5 -2 12.2168303 -0.36334223 1.08110221
12 -5 -1 12.1421622 -0.22454581 1.43526214
13 -4 -1 10.5502143 -0.86015884 0.62132264
14 -3 -1 9.5431203 -0.55845539 0.58204048
15 -2 -1 8.6135278 0.58556445 1.03000890
16 -1 -1 5.9954109 1.79818632 1.13151965
17 0 -1 2.2277653 0.84147098 0.84147098
18 1 -1 5.0861135 -0.88888889 0.22222222
19 2 -1 7.6309368 -1.28591543 0.04741790
20 3 -1 9.3906254 -1.10764453 0.42954555
21 4 -1 11.3680186 -0.23988544 1.43912691
22 5 -1 13.1783802 0.70193281 2.47148014
23 -5 0 11.4744286 -0.95892427 0.96415265
24 -4 0 10.7781363 -0.75680250 1.12555045
25 -3 0 10.5062902 0.14112001 1.94112001
26 -2 0 9.5681733 0.90929743 2.50929743
27 -1 0 6.2277653 0.84147098 1.84147098
28 0 0 0.0000000 0.00000000 0.00000000 ### check
29 1 0 4.5448234 -0.84147098 0.15852902
30 2 0 7.7495784 -0.90929743 0.69070257
31 3 0 10.2240502 -0.14112001 1.65887999
32 4 0 12.2917413 0.75680250 2.63915544
33 5 0 13.3922771 0.95892427 2.88200120
34 -5 1 11.6647752 -0.81167218 0.95787515
35 -4 1 11.6502586 0.04235458 1.72136692
36 -3 1 11.2092202 0.71095032 2.24814040
37 -2 1 9.3138788 0.39702654 1.73035987
38 -1 1 5.0861135 -0.88888889 0.22222222
39 0 1 0.5448234 -0.84147098 -0.84147098
40 1 1 4.1768160 -0.02040854 -0.68707520
41 2 1 8.3312878 0.30332444 0.74776888
42 3 1 11.0567253 0.95514960 2.09564547
43 4 1 12.4680629 1.05768971 2.53917119
44 5 1 12.7009932 0.33428518 1.99409314
45 -5 2 12.4990703 -0.08110221 1.36334223
46 -4 2 12.2953945 0.54648564 1.85260808
47 -3 2 10.8950958 0.22922609 1.39249139
48 -2 2 8.0487701 -0.98765432 0.09876543
49 -1 2 4.2976035 -1.95258210 -0.95258210
50 0 2 2.3095784 -0.90929743 -0.90929743
51 1 2 4.9979545 0.96999110 -0.25223112
52 2 2 8.8055726 1.74445682 0.85556793
53 3 2 11.0125491 1.57116917 1.50994468
54 4 2 11.6655125 0.64222729 1.22272616
55 5 2 11.7009637 -0.43476438 0.56523562
>
Bill Venables
CSIRO Laboratories
PO Box 120, Cleveland, 4163
AUSTRALIA
Office Phone (email preferred): +61 7 3826 7251
Fax (if absolutely necessary): +61 7 3826 7304
Mobile: +61 4 8819 4402
Home Phone: +61 7 3286 7700
mailto:[EMAIL PROTECTED]
http://www.cmis.csiro.au/bill.venables/
-----Original Message-----
From: [EMAIL PROTECTED] [mailto:[EMAIL PROTECTED]
On Behalf Of Lavan
Sent: Tuesday, 25 March 2008 10:00 AM
To: r-help@r-project.org
Subject: [R] derivatives in R
Hi, I posted this message earlier in "Rmetrics" and I don't know whether
I
posted in the wrong place, so I'm posting it again in Rhelp.
I have a function in x and y and let's call it f(x,y). I need to get the
Hessian matrix. i.e I need (d^2f/dx^2), (d^2f/dxdy), (d^2f/dydx),
(d^2f/dy^2). I can get these using the D function. now I need to
evaluste the
hessian matrix for -5<x<5 and -2<y<2 and I also need to print the
output.
I tried using a for loop, but it is giving me an error.
Thanks
--
View this message in context:
http://www.nabble.com/derivatives-in-R-tp16265419p16265419.html
Sent from the R help mailing list archive at Nabble.com.
______________________________________________
R-help@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@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.
