Hello R users,
I am trying to estimate the parameters of a bimodal normal distribution using
moments matching, so I have to solve a non-linear system of equations. How can I
solve the following simple example?
x^2 - y^2 = 6
x ? y = 3
I heard about nlsystemfit, but I don?t know how to run it exactly. I have tried
the following code, but it doesn?t really work:
f1 <-y~ x[1]^2-x[2]^2-6
f2 <-z~ x[1]-x[2]-3
f <- list(f1=0,f2=0)
nlsystemfit("OLS",f,startvals=c(0,0))
Thank You in advance for your help.
Evgeniq Petrova
2008/4/25 Radka Pancheva <radica at abv.bg>:> I am trying to estimate the parameters of a bimodal normal distribution using moments matching, so I have to solve a non-linear system of equations. How can I solve the following simple example? > > x^2 - y^2 = 6 > x ? y = 3 > > I heard about nlsystemfit, but I don't know how to run it exactly. I have tried the following code, but it doesn't really work: > > > f1 <-y~ x[1]^2-x[2]^2-6 > f2 <-z~ x[1]-x[2]-3 > f <- list(f1=0,f2=0) > nlsystemfit("OLS",f,startvals=c(0,0))You could try the recent package BB by Ravi Varadhan. The code could be the following: library(BB) f <- function(x) { x1 <- x[1] x2 <- x[2] F <- rep(NA, 2) F[1] <- x1^2 - x2^2 - 6 F[2] <- x1 - x2 - 3 return(F) } p0 <- c(1,2) dfsane(par=p0, fn=f,control=list(maxit=3000)) I got the solution: x1 = 2.5 x2 = -0.5 Paul
Hello Paul,
Thank you for your quick answer. I have tried to use your advice and to estimate
the parameters of beta distribution with moments matching. This is my code:
ex <- 0.3914877
ex2 <- 0.2671597
my.mm <- function(x){
p <- x[1]
q <- x[2]
p <- .Machine$double.eps
q <- .Machine$double.eps
F <- rep(NA,2)
F[1] <- p/(p + q)
F[2]<- (p*q + (p + q + 1)*p^2)/((p + q + 1)*(p + q)^2)
return(F)
}
p0 <- c(ex,ex2)
dfsane(par=p0, fn=my.mm,control=list(maxit=50000))
and I became the following output:
?
iteration: 3640 ||F(xn)|| = 0.7071068
iteration: 3641 ||F(xn)|| = 0. 7071068
?
iteration: 49990 ||F(xn)|| = 0. 7071068
iteration: 50000 ||F(xn)|| = 0. 7071068
$par
[1] -446.2791 -446.4034
$residual
[1] 0.5
$fn.reduction
[1] 0
$feval
[1] 828495
$iter
[1] 50001
$convergence
[1] 1
$message
[1] "Maximum limit for iterations exceeded"
I have tried maxiter=100000 but the output is the same. I know that ex and ex2
are bringing the problems, but I am stuck with them. How can I make it
convergent?
Thank you,
Evgeniq
>2008/4/25 Radka Pancheva <radica at abv.bg>:
>> I am trying to estimate the parameters of a bimodal normal
distribution using moments matching, so I have to solve a non-linear system of
equations. How can I solve the following simple example?
>>
>> x^2 - y^2 = 6
>> x ? y = 3
>>
>> I heard about nlsystemfit, but I don't know how to run it
exactly. I have tried the following code, but it doesn't really work:
>>
>>
>> f1 <-y~ x[1]^2-x[2]^2-6
>> f2 <-z~ x[1]-x[2]-3
>> f <- list(f1=0,f2=0)
>> nlsystemfit("OLS",f,startvals=c(0,0))
>
>You could try the recent package BB by Ravi Varadhan. The code could
>be the following:
>
>library(BB)
>
>f <- function(x) {
> x1 <- x[1]
> x2 <- x[2]
>
> F <- rep(NA, 2)
>
> F[1] <- x1^2 - x2^2 - 6
> F[2] <- x1 - x2 - 3
>
> return(F)
>}
>
>p0 <- c(1,2)
>dfsane(par=p0, fn=f,control=list(maxit=3000))
>
>I got the solution:
>
>x1 = 2.5
>x2 = -0.5
>
>Paul
>
>______________________________________________
>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.
>
> x^2 - y^2 = 6 > x ? y = 3 > >You can also try this # function f <- function(x) { y <- numeric(2) y[1] <- x[1]^2-x[2]^2-6 y[2] <- x[1]-x[2]-3 y } # function values transformed to scalar # minimising fnorm this way is not the best method of finding a solution for f(x)=0 # there may be values for x which minimise fnorm but do not set f(x) = 0 # but you can always try fnorm <- function(z) {p <- f(z);return(crossprod(p))} #starting values xstart <- c(0,0) # You can use nlm or nlminb nlm(fnorm,xstart) nlminb(xstart,fnorm) Sometimes minpack.lm can be used to find roots. Do library(minpack.lm) nls.lm(xstart,f) and in this case you'll see that it doesn't work. In this case stick to nlm and/or nlminb. Berend -- View this message in context: http://www.nabble.com/Non-linear-system-of-equations-tp16893056p16895541.html Sent from the R help mailing list archive at Nabble.com.
Try this:> library(Ryacas) > x <- Sym("x") > y <- Sym("y") > Solve(List(x^2+y^2==6, x-y==3), List(x,y))expression(list(list(x == root(6 - y^2, 2), y == y))) 2008/4/25 Radka Pancheva <radica at abv.bg>:> Hello R users, > > I am trying to estimate the parameters of a bimodal normal distribution using moments matching, so I have to solve a non-linear system of equations. How can I solve the following simple example? > > x^2 - y^2 = 6 > x ? y = 3 > > I heard about nlsystemfit, but I don't know how to run it exactly. I have tried the following code, but it doesn't really work: > > > f1 <-y~ x[1]^2-x[2]^2-6 > f2 <-z~ x[1]-x[2]-3 > f <- list(f1=0,f2=0) > nlsystemfit("OLS",f,startvals=c(0,0)) > > Thank You in advance for your help. > > Evgeniq Petrova > > ______________________________________________ > 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. >
Hello,
Thank you all for your answers. I'll try them all and I'm sure it will
work.
>
>Try this:
>
>> library(Ryacas)
>> x <- Sym("x")
>> y <- Sym("y")
>> Solve(List(x^2+y^2==6, x-y==3), List(x,y))
>expression(list(list(x == root(6 - y^2, 2), y == y)))
>
>
>2008/4/25 Radka Pancheva <radica at abv.bg>:
>> Hello R users,
>>
>> I am trying to estimate the parameters of a bimodal normal
distribution using moments matching, so I have to solve a non-linear system of
equations. How can I solve the following simple example?
>>
>> x^2 - y^2 = 6
>> x ? y = 3
>>
>> I heard about nlsystemfit, but I don't know how to run it exactly.
I have tried the following code, but it doesn't really work:
>>
>>
>> f1 <-y~ x[1]^2-x[2]^2-6
>> f2 <-z~ x[1]-x[2]-3
>> f <- list(f1=0,f2=0)
>> nlsystemfit("OLS",f,startvals=c(0,0))
>>
>> Thank You in advance for your help.
>>
>>
>> ______________________________________________
>> 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.
>>
>