hi thanks,
pls see the code below. the elements of matrix Il should be the derivatives
in matrix h (in Boldface) plus some additional terms, here I did calculate
derivatoves in the matrix, but I want R to calculate it for me.
pid<- function (tau1,tau2,sigma=10,n=100,x_bar=20) {
library(akima)
library(lattice)
library(rgl)
a1 <- 1.5
a2 <- 4
t1 <-seq(-1,1);
t2 <-seq(3,4);
l <-length(t1);
m <-length(t2);
N1 <- matrix(NA,l,m);
N2 <- matrix(NA,l,m);
wp <- matrix(NA,l,m);
rr <- matrix(NA,l,m);
N <- matrix(NA,l,m);
mu <- matrix(NA,l,m);
Il <- matrix(NA,l,m);
iIl <- matrix(NA,l,m);
theta <-expression(x+y);
dx <-D(theta, "x");
dy <-D(theta, "y");
d2x <-D(dx, "x");
dxy <-D(dx, "y");
dyx <-D(dy, "x");
d2y <-D(dy, "y");
h <- matrix(c(d2x,dxy,dyx,d2y),nrow=2,ncol=2);
for (i in 1:l)
{
for (j in 1:m)
{
wp[i,j]<-((t1[i]-a1)/tau1)^2+((t2[j]-a2)/tau2)^2
mu[i,j] <-(t1[i] + t2[j]);
Il[i,j]
<-matrix(c(n/(sigma^2)+1/(tau1^2),n/(sigma^2),n/(sigma^2),n/(sigma^2)+1/
tau2^2)),2,2);
iIl <-solve(Il);
ul[i,j] <-
(c(n*(x_bar-mu[i,j])/(sigma^2)-(t1[i]-a1)/(tau1^2),n*(x_bar-mu[i,j])/
(sigma^2)-(t2[j]-a2)/(tau2^2)));
ul <-matrix(c(n*(x_bar-mu)/(sigma^2)-(theta1-a1)/(tau1^2),n*(x_bar-mu)/
(sigma^2)-(theta2-a2)/(tau2^2)),2,1);
ult <-t(ul);
wl[i,j] <-t(ul)%*%solve(Il,ul);
rr[i,j] <-(wp[i,j]/wl[i,j]);
N[i,j] <-(wl[i,j]-wp[i,j])/wl[i,j]
}
}
list("Test Statistics Wp"= wp,"Stat Ratio"=rr, "N"=N,"N1"=N1,
"N2"=N2, "mean"=mu, "Il"=Il,
"inverse of Il"=iIl)
}
--
View this message in context:
http://www.nabble.com/derivatives-in-R-tp16265419p16267397.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.
