[R] derivatives in R

Tue Mar 25 04:00:07 CET 2008

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:Bill.Venables at csiro.au
http://www.cmis.csiro.au/bill.venables/ 

-----Original Message-----
From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org]
On Behalf Of Lavan
Sent: Tuesday, 25 March 2008 10:00 AM
To: r-help at 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 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.