search for: jacobian

...arameters to estimate      alp <- pars[1]    b1  <- pars[2]     b2  <- pars[3]    rho <- pars[4]   f1 <- pars[1]*pars[2] f2 <- pars[1]*(pars[1]+1)*pars[2]^2 f3 <- pars[1]*(pars[1]+1)*pars[3]^2 f4 <- pars[1]*pars[2]*pars[3](pars[1]+pars[4]) fval  <- c(f1,f2,f3,f4) ## jacobian matrix j11 <- pars[2] j12 <- pars[1] j13 <- 0 j14 <- 0 j21 <-  (2*pars[1]+1)*pars[2]^2 j22 <-  2*pars[2]*pars[1]*(pars[1]+1) j23 <-  0 j24 <-  0 j31 <- (2*pars[1]+1)*pars[3]^2 j32 <- 2*pars[3]*pars[1]*(pars[1]+1) j33 <- 0 j34 <- 0 j41 <- pars[2]*pars[3]*(2*pa...

Hi, I am trying to solve a model that consists of rather stiff ODEs in R. I use the package ODEsolve (lsoda) to solve these ODEs. To speed up the integration, the jacobian is also specified. Basically, the model is a one-dimensional advection-diffusion problem, and thus the jacobian is a tridiagonal matrix. The size of this jacobian is 100*100. In the original package LSODA it is possible to specify that the jacobian is banded, which makes its inversion ve...

...alling C code. > > > > Any thoughts on how to best handle the situation? Is this a bug in dpoco? Is there a simple way to test for any NaNs in a vector? > You should/could use macro R_FINITE to test each entry of the hessian. > In package nleqslv I test for a "correct" jacobian like this in file nleqslv.c in function fcnjac: > for (j = 0; j < *n; j++) > for (i = 0; i < *n; i++) { > if( !R_FINITE(REAL(sexp_fjac)[(*n)*j + i]) ) > error("non-finite value(s) returned by jacobian (row=%d,col=%d)",i+1,j+1); &g...

...t;>> Any thoughts on how to best handle the situation? Is this a bug in dpoco? Is there a simple way to test for any NaNs in a vector? >> >>> You should/could use macro R_FINITE to test each entry of the hessian. >>> In package nleqslv I test for a "correct" jacobian like this in file nleqslv.c in function fcnjac: >> >>> for (j = 0; j < *n; j++) >>> for (i = 0; i < *n; i++) { >>> if( !R_FINITE(REAL(sexp_fjac)[(*n)*j + i]) ) >>> error("non-finite value(s) returned by jacobi...

...th nls.lm function of minpack.lm package. For the first fit, the parameter estimates keep changing even after 1000 iterations (Th) and I have a following error message for fit of hydraulic conductivity (k); Reason for termination: The cosine of the angle between `fvec' and any column of the Jacobian is at most `gtol' in absolute value. I dont know what I can do for solve my problem Yesterday I asked to Katharine Mullen Here are her aswer Yann P?riard ?tudiant ? la ma?trise en sols et environnement D?partement des sols et de g?nie agroalimentaire Centre de Recherche en Horticulture...

In my package 'glmmML' I'm using old C code and linpack in the optimizing procedure. Specifically, one part of the code looks like this: F77_CALL(dpoco)(*hessian, &bdim, &bdim, &rcond, work, info); if (*info == 0){ F77_CALL(dpodi)(*hessian, &bdim, &bdim, det, &job); ........ This usually works OK, but with an ill-conditioned data

...rentiation: library(mvtnorm) library(numDeriv) f=function(pars, xx, yy) { mu=pars[1:2] sig1=pars[3] sig2=pars[4] rho=pars[5] sig=matrix(c(sig1^2,rho*sig1*sig2,rho*sig1*sig2,sig2^2),2) dmvnorm(cbind(x,y),mu,sig) } mu1=1 mu2=2 sig1=3 sig2=4 rho=.5 x=2 # or a x vector y=3 # or a y vector jacobian(f,c(mu1,mu2,sig1,sig2,rho),xx=x,yy=y) # (Can replace ?jacobian? with ?grad? if x and y have length 1.) -- Karl Ove Hufthammer

I'm trying to obtain numerical derivative of a probability computed with mvtnorm with respect to its parameters using grad() and jacobian() from NumDeriv. To simplify the matter, here is an example: PP1 <- function(p){ thetac <- p thetae <- 0.323340333 thetab <- -0.280970036 thetao <- 0.770768082 ssigma <- diag(4) ssigma[1,2] <- 0.229502120 ssigma[1,3] <- 0.677949335 ssigma[1,4] &...

...00 386045.812 379593.1 39762.40 [4,] 420.000 42007.764 39762.4 5740.00 $code [1] 1 $iterations [1] 31 I was under the impression that you could also obtain the se of the parameter estimates using the sqrt( diag( res$hessian ) ), but I haven't been able to reproduce the se computed by the Jacobian se <- sqrt( mse * diag( solve( crossprod( J ) ) ) ) # gives the correct results... hse <- sqrt( ( res$minimum / 8 ) * diag( solve( res$hessian ) ) ) # gives similar results, but why 8? I've tried to put the functionality to include the jacobian and hessian in the objective function...

Hi R-users, Does R has a topic on newton's method? Thank you for the info.

...t;>> >>> Any thoughts on how to best handle the situation? Is this a bug in dpoco? Is there a simple way to test for any NaNs in a vector? > >> You should/could use macro R_FINITE to test each entry of the hessian. >> In package nleqslv I test for a "correct" jacobian like this in file nleqslv.c in function fcnjac: > >> for (j = 0; j < *n; j++) >> for (i = 0; i < *n; i++) { >> if( !R_FINITE(REAL(sexp_fjac)[(*n)*j + i]) ) >> error("non-finite value(s) returned by jacobian (row=%d,col=%d...

I've been trying to figure out how to build a list of terms from a nonlinear model (terms() returns a error). I need to compute and evaluate the partial derivatives (Jacobian) for each equaiton in a set of equations. For example: > eqn <- q ~ s0 + s1 * p + s2 * f + s3 * a > sv2 <- c(d0=3,d1=4.234,d2=4,s0=-2.123,s1=0.234,s2=2.123,s3=4.234) > names( sv2 ) [1] "d0" "d2" "d1" "s0" "s2" "s3" "...

I am using nonlinear regression to fit a couple of variables to a set of measurements. I would like to do some significance tests for the estimated parameters. I am able to check the confidence intervals using the Jacobian coming out of nonlinear regression. I do see in a paper which shows t-value (it says estimated by White method??), f-value, f-test, and j-test, are these available in matlab, or I could code them myself but I need details about these tests. I would be glad if you could provide a reference... Tha...

Hi R-users, Can Maple function be exported to R? I have a jacobian matrix (4X4)?from maple in algebraic form which involve modified Bessel function of the first kind. I just wonder whether we can use algebraic form into R before the value of the parameters can be estimated. Thank you so much for your attention and help.

Dear friends - I have a function for the charge in a fluid (water) buffered with HEPES and otherwise only containing Na and Cl so that [Na] - [Cl] = SID (strong ion difference) goes from -1 mM to 1 mM. With known SID and total HEPES concentration I can calculate accurately the pH if I know 3 pK values for HEPES by finding the single root with uniroot Now, the problem is that there is some

Hello The excellent adapt package integrates over multi-dimensional hypercubes. I want to integrate over a multidimensional simplex. Has anyone implemented such a thing in R? I can transform an n-simplex to a hyperrectangle but the Jacobian is a rapidly-varying (and very lopsided) function and this is making adapt() slow. [ A \dfn{simplex} is an n-dimensional analogue of a triangle or tetrahedron. It is the convex hull of (n+1) points in an n-dimensional Euclidean space. My application is a variant of the Dirichlet distribution:...

Dear All, I couple of weeks ago, I’ve asked for a package recommendation for nonlinear optimization. In my problem I have a fairly complicated non-linear objective function subject to one non-linear equality constrain. I’ve been suggested to use the *Rdonlp2* package, but I did not get any results after running the program for 5 hrs. Is it normal to run this type of programs for hours? Also,

...problem in which first and second derivatives are computed numerically. This is part of a numerical method course I am teaching in which I want to compare speeds of R versus Fortran (We solve a general equilibrium problem all numerically, if you want to know). Because of this discrepancy, the Jacobian and Hessian in R versus in Fortran are quite different, which results in the Newton method producing a different solution (for a given stopping rule). Since the solution computed in Fortran is almost identical to the analytical solution, I suspect that the sum in Fortran may be more accurate (T...

...t;\nWarning: weights matrix not in sparse form\n") switch(type, lag = if (!quiet) cat("\nSpatial lag model\n"), mixed = if (!quiet) cat("\nSpatial mixed autoregressive model\n"), stop("\nUnknown model type\n")) if (!quiet) cat("Jacobian calculated using weights matrix eigenvalues\n") y <- model.response(mf, "numeric") if (any(is.na(y))) stop("NAs in dependent variable") x <- model.matrix(mt, mf) if (any(is.na(x))) stop("NAs in independent variable") if (nrow(x) != nrow(w...

SSweibull() :  problems with step factor and singular gradient Hello I am working with growth data of ~4000 tree seedlings and trying to fit non-linear Weibull growth curves through the data of each plant. Since they differ a lot in their shape, initial parameters cannot be set for all plants. That’s why I use the self-starting function SSweibull(). However, I often got two error messages: