[R] Solve system of non linear equations using nasted loops

[Berend Hasselman]

Your code looks needlessly complicated and inefficient in my view.
Even if you don't need or use nested for loops it can still be simplified.

1: First define  some new constants to avoid the  the M? stuff in your function

<code>
pME <- ME[2:9]/ME[1]
pMI <- MI[2:9]/MI[1]
pMA <- MA[2:9]/MA[1]
pMAE <- (MA[1]-ME[1])/(MA[2:9]-ME[2:9])
pMAI <- (MA[2:9]-MI[2:9])/(MA[1]-MI[1])
pMIA <- (MI[2:9]-MA[2:9])/(MI[1]-MA[1])
</code>

2. Redefine your function and avoid silly and unnecessary  use of 1* and 2* and use a constant for log(.../...)
The vector pXX is a vector consisting of elements from the vector pME etc.

<code>
fun <- function(x,pXX) {
    f <- numeric(length(x)) # read as:
    pME <- pXX[1]
    pMI <- pXX[2]
    pMA <- pXX[3]
    pMAE <- pXX[4]
    pMAI <- pXX[5]
    pMIA <- pXX[6]
    const <- log(0.19245009/0.384900179)
    f[1] <- x[3] - (log(pME)+(x[1]+x[2])*2*(log(pMAE)))/const
    f[2] <- x[3] - (log(pMI)-(x[2])*2*(log(pMAI)))/const
    f[3] <- x[3] - (log(pMA)-(x[1])*2*(log(pMIA)))/const
    f
}
</code>

This can be simplified further by doing the log(pXX) beforehand. 

3. I'm not going to have a go at using expand.grid. Too lazy and more involved to get everything correct.
So use a nested loop as follows:

- prepare a result matrix of correct dimensions with7 columns: one for nleaslv termination code, 3 for the result x-vales and another 3 to store i,j,k.
- storing the termination code is a must to see if convergence was obtained.
- and before invoking nleqslv test if fun(startx) contains NA values. If it does fill the appropriate row of Result with NA.
- use an index counter (indx) to store the results in the rows of Result

<code>
Result <- matrix(0,nrow=9*9*9,ncol=7)

indx <- 1
for (i in 1:9) {
    for (j in 1:9) {
        for (k in 1:9) {
            pXXijk <- c(pME[i],pMI[j],pMA[k],pMAE[i],pMAI[j],pMIA[k])
            f.startx <- fun(startx,pXXijk)
            if(anyNA(f.startx)) {
                Result[indx,1:4] <- NA
            } else {
                z <- nleqslv(startx,fun,pXX=pXXijk)
                Result[indx,1:4] <- c(z$termcd,z$x)
            }
            Result[indx,5:7] <- c(i,j,k)
            indx <- indx+1
        }
    }
}
</code>

And now you can inspect Result, which contains 729 rows (and not 1000) and continue to make  a data.frame with column names.
The rest is up to you.

Berend

On 10 May 2017, at 20:20, Santiago Burone <santiago
burone at icloud.com> wrote:
> 
> Hello, 
> 
> I'm new at R and I would like to use it in order to solve a system of non linear equations. I have the code that works but im not able to save the results. 
> 
> 
> 
> My system has three equations and i would like to solve this using three nested loops, becouse i need solutions for all the values combinations.
> 
> 
> 
> 
> 
> The code im using is this:
> 
> 
> 
> library(nleqslv)
> #x3 es gamma
> #x2 es alpha
> #x1 es beta
> MA <- c(50000, 43600, 40000, 38800, 37600, 34400, 31600, 27200, 24400, 20000)
> MI <- c(10000, 21800, 20000, 19400, 18800, 17200, 15800, 13600, 12200, 10000)
> ME <- c(30000, 32700, 30000, 29100, 28200, 25800, 23700, 20400, 18300, 15000)
> DE <- c(0.384900179, 0.19245009, 0.19245009, 0.19245009, 0.19245009, 0.19245009, 0.19245009, 0.19245009, 0.19245009, 0.19245009)
> for (i in 1:9) {
> for (j in 1:9){
> for (k in 1:9){
> df.names <- paste("SociedadB",1:10,sep="")
> 
> fun <- function(x) { 
> f <- numeric(length(x)) # read as:
> f[1] <- 1*x[3] - (log(ME[1+i]/ME[1])+(1*x[1]+1*x[2])*(2)*(log((MA[1]-ME[1])/(MA[i+1]-ME[i+1]))))/(log(0.19245009/0.384900179))
> f[2] <- 1*x[3] - (log(MI[1+j]/MI[1])-(1*x[2])*(2)*(log((MA[1+j]-MI[1+j])/(MA[1]-MI[1]))))/(log(0.19245009/0.384900179)) 
> f[3] <- 1*x[3] - (log(MA[1+k]/MA[1])-(1*x[1])*(2)*(log((MI[1+k]-MA[1+k])/(MI[1]-MA[1]))))/(log(0.19245009/0.384900179)) 
> f 
> } 
> startx <- c(1,1,1) 
> answers<-as.data.frame(nleqslv(startx,fun)) 
> d.frame <- (answers$x)
> assign(df.names[k+20], d.frame)
> despejes<- data.frame(SociedadB2, SociedadB3, SociedadB4, SociedadB5, SociedadB6, SociedadB7, SociedadB8, SociedadB9, SociedadB10)
> 
> }
> }
> }
> 
> despejes<- data.frame(SociedadB2, SociedadB3, SociedadB4, SociedadB5, SociedadB6, SociedadB7, SociedadB8, SociedadB9, SociedadB10)
> 
> 
> 
> If im not wrong i should have 1000 different solutions but the data frame only save 10 results. I tried copying the last part of the code between the last two brackets (changing the name of the data frame despejes for despejes 1, 2 and 3 but i only get three data frames with exactly the same results. I dont know how could i save all the results of every combination. I mean, i need to solve for i=1, j=1 and the 9 diferent values of k and save those results, then i=1, j=2 and the 9 values of k..
> 
> 
> 
> I know this might be a easy question but i have not been able to find a solution in the last days. 
> 
> 
> 
> Thanks in advance if you can help me.
> 
> 
> Santiago. 
> 	[[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.