Re: [R] solving integral equations with undefined parameters using multiroot
lm (cl$y ~ cl$x)$coef (Intercept)cl$x 0.1817509 -1.000 On Fri, May 7, 2021 at 1:56 PM Abbs Spurdle wrote: > > #using vF1 function > #from my previous posts > u <- seq (0, 0.25,, 200) > cl <- contourLines (u, u, outer (u, u, vF1),, 0)[[1]] > plot (cl$x, cl$y, type="l") > > > On Thu, May 6, 2021 at 10:18 PM Ursula Trigos-Raczkowski > wrote: > > > > Thanks for your reply. Unfortunately the code doesn't work even when I > > change the parameters to ensure I have "different" equations. > > Using mathematica I do see that my two equations form planes, intersecting > > in a line of infinite solutions but it is not very accurate, I was hoping R > > would be more accurate and tell me what this line is, or at least a set of > > solutions. > > > > On Thu, May 6, 2021 at 5:28 AM Abbs Spurdle wrote: > >> > >> Just realized five minutes after posting that I misinterpreted your > >> question, slightly. > >> However, after comparing the solution sets for *both* equations, I > >> can't see any obvious difference between the two. > >> If there is any difference, presumably that difference is extremely small. > >> > >> > >> On Thu, May 6, 2021 at 8:39 PM Abbs Spurdle wrote: > >> > > >> > Hi Ursula, > >> > > >> > If I'm not mistaken, there's an infinite number of solutions, which > >> > form a straight (or near straight) line. > >> > Refer to the following code, and attached plot. > >> > > >> > begin code--- > >> > library (barsurf) > >> > vF1 <- function (u, v) > >> > { n <- length (u) > >> > k <- numeric (n) > >> > for (i in seq_len (n) ) > >> > k [i] <- intfun1 (c (u [i], v [i]) ) > >> > k > >> > } > >> > plotf_cfield (vF1, c (0, 0.2), fb = (-2:2) / 10, > >> > main="(integral_1 - 1)", > >> > xlab="S[1]", ylab="S[2]", > >> > n=40, raster=TRUE, theme="heat", contour.labels=TRUE) > >> > end code > >> > > >> > I'm not familiar with the RootSolve package. > >> > Nor am I quite sure what you're trying to compute, given the apparent > >> > infinite set of solutions. > >> > > >> > So, for now at least, I'll leave comments on the root finding to someone > >> > who is. > >> > > >> > > >> > Abby > >> > > >> > > >> > On Thu, May 6, 2021 at 8:46 AM Ursula Trigos-Raczkowski > >> > wrote: > >> > > > >> > > Hello, > >> > > I am trying to solve a system of integral equations using multiroot. I > >> > > have > >> > > tried asking on stack exchange and reddit without any luck. > >> > > Multiroot uses the library(RootSolve). > >> > > > >> > > I have two integral equations involving constants S[1] and S[2] (which > >> > > are > >> > > free.) I would like to find what *positive* values of S[1] and S[2] > >> > > make > >> > > the resulting > >> > > (Integrals-1) = 0. > >> > > (I know that the way I have the parameters set up the equations are > >> > > very > >> > > similar but I am interested in changing the parameters once I have the > >> > > code > >> > > working.) > >> > > My attempt at code: > >> > > > >> > > ```{r} > >> > > a11 <- 1 #alpha_{11} > >> > > a12 <- 1 #alpha_{12} > >> > > a21 <- 1 #alpha_{21} > >> > > a22 <- 1 #alpha_{22} > >> > > b1 <- 2 #beta1 > >> > > b2 <- 2 #beta2 > >> > > d1 <- 1 #delta1 > >> > > d2 <- 1 #delta2 > >> > > g <- 0.5 #gamma > >> > > > >> > > > >> > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > >> > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > >> > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > >> > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > >> > > > >> > > #defining equation we would like to solve > >> > > intfun1<- function(S) {integrate(function(x) integrand1(x, > >> > > S),lower=0,upper=Inf)[[1]]-1} > >> > > intfun2<- function(S) {integrate(function(x) integrand2(x, > >> > > S),lower=0,upper=Inf)[[1]]-1} > >> > > > >> > > #putting both equations into one term > >> > > model <- function(S) c(F1 = intfun1,F2 = intfun2) > >> > > > >> > > #Solving for roots > >> > > (ss <-multiroot(f=model, start=c(0,0))) > >> > > ``` > >> > > > >> > > This gives me the error Error in stode(y, times, func, parms = parms, > >> > > ...) : > >> > > REAL() can only be applied to a 'numeric', not a 'list' > >> > > > >> > > However this simpler example works fine: > >> > > > >> > > ```{r} > >> > > #Defining the functions > >> > > model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) > >> > > > >> > > #Solving for the roots > >> > > (ss <- multiroot(f = model, start = c(0,0))) > >> > > ``` > >> > > > >> > > Giving me the required x_1= 4 and x_2 =1. > >> > > > >> > > I was given some code to perform a least squares analysis on the same > >> > > system but I neither understand the code, nor believe that it is doing > >> > > what > >> > > I am looking for as different initial values give wildly different S > >> > > values. > >> > > > >> > > ```{r} > >> > > a11 <- 1 #alpha_{11} > >> > > a12 <- 1 #alpha_{12} > >> > > a21 <- 1 #alpha_{21} > >> > > a22 <- 1
Re: [R] solving integral equations with undefined parameters using multiroot
#using vF1 function #from my previous posts u <- seq (0, 0.25,, 200) cl <- contourLines (u, u, outer (u, u, vF1),, 0)[[1]] plot (cl$x, cl$y, type="l") On Thu, May 6, 2021 at 10:18 PM Ursula Trigos-Raczkowski wrote: > > Thanks for your reply. Unfortunately the code doesn't work even when I change > the parameters to ensure I have "different" equations. > Using mathematica I do see that my two equations form planes, intersecting in > a line of infinite solutions but it is not very accurate, I was hoping R > would be more accurate and tell me what this line is, or at least a set of > solutions. > > On Thu, May 6, 2021 at 5:28 AM Abbs Spurdle wrote: >> >> Just realized five minutes after posting that I misinterpreted your >> question, slightly. >> However, after comparing the solution sets for *both* equations, I >> can't see any obvious difference between the two. >> If there is any difference, presumably that difference is extremely small. >> >> >> On Thu, May 6, 2021 at 8:39 PM Abbs Spurdle wrote: >> > >> > Hi Ursula, >> > >> > If I'm not mistaken, there's an infinite number of solutions, which >> > form a straight (or near straight) line. >> > Refer to the following code, and attached plot. >> > >> > begin code--- >> > library (barsurf) >> > vF1 <- function (u, v) >> > { n <- length (u) >> > k <- numeric (n) >> > for (i in seq_len (n) ) >> > k [i] <- intfun1 (c (u [i], v [i]) ) >> > k >> > } >> > plotf_cfield (vF1, c (0, 0.2), fb = (-2:2) / 10, >> > main="(integral_1 - 1)", >> > xlab="S[1]", ylab="S[2]", >> > n=40, raster=TRUE, theme="heat", contour.labels=TRUE) >> > end code >> > >> > I'm not familiar with the RootSolve package. >> > Nor am I quite sure what you're trying to compute, given the apparent >> > infinite set of solutions. >> > >> > So, for now at least, I'll leave comments on the root finding to someone >> > who is. >> > >> > >> > Abby >> > >> > >> > On Thu, May 6, 2021 at 8:46 AM Ursula Trigos-Raczkowski >> > wrote: >> > > >> > > Hello, >> > > I am trying to solve a system of integral equations using multiroot. I >> > > have >> > > tried asking on stack exchange and reddit without any luck. >> > > Multiroot uses the library(RootSolve). >> > > >> > > I have two integral equations involving constants S[1] and S[2] (which >> > > are >> > > free.) I would like to find what *positive* values of S[1] and S[2] make >> > > the resulting >> > > (Integrals-1) = 0. >> > > (I know that the way I have the parameters set up the equations are very >> > > similar but I am interested in changing the parameters once I have the >> > > code >> > > working.) >> > > My attempt at code: >> > > >> > > ```{r} >> > > a11 <- 1 #alpha_{11} >> > > a12 <- 1 #alpha_{12} >> > > a21 <- 1 #alpha_{21} >> > > a22 <- 1 #alpha_{22} >> > > b1 <- 2 #beta1 >> > > b2 <- 2 #beta2 >> > > d1 <- 1 #delta1 >> > > d2 <- 1 #delta2 >> > > g <- 0.5 #gamma >> > > >> > > >> > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* >> > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} >> > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* >> > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} >> > > >> > > #defining equation we would like to solve >> > > intfun1<- function(S) {integrate(function(x) integrand1(x, >> > > S),lower=0,upper=Inf)[[1]]-1} >> > > intfun2<- function(S) {integrate(function(x) integrand2(x, >> > > S),lower=0,upper=Inf)[[1]]-1} >> > > >> > > #putting both equations into one term >> > > model <- function(S) c(F1 = intfun1,F2 = intfun2) >> > > >> > > #Solving for roots >> > > (ss <-multiroot(f=model, start=c(0,0))) >> > > ``` >> > > >> > > This gives me the error Error in stode(y, times, func, parms = parms, >> > > ...) : >> > > REAL() can only be applied to a 'numeric', not a 'list' >> > > >> > > However this simpler example works fine: >> > > >> > > ```{r} >> > > #Defining the functions >> > > model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) >> > > >> > > #Solving for the roots >> > > (ss <- multiroot(f = model, start = c(0,0))) >> > > ``` >> > > >> > > Giving me the required x_1= 4 and x_2 =1. >> > > >> > > I was given some code to perform a least squares analysis on the same >> > > system but I neither understand the code, nor believe that it is doing >> > > what >> > > I am looking for as different initial values give wildly different S >> > > values. >> > > >> > > ```{r} >> > > a11 <- 1 #alpha_{11} >> > > a12 <- 1 #alpha_{12} >> > > a21 <- 1 #alpha_{21} >> > > a22 <- 1 #alpha_{22} >> > > b1 <- 2 #beta1 >> > > b2 <- 2 #beta2 >> > > d1 <- 1 #delta1 >> > > d2 <- 1 #delta2 >> > > g <- 0.5 #gamma >> > > >> > > >> > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* >> > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} >> > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* >> > >
Re: [R] solving integral equations with undefined parameters using multiroot
Thanks for your reply. Unfortunately the code doesn't work even when I change the parameters to ensure I have "different" equations. Using mathematica I do see that my two equations form planes, intersecting in a line of infinite solutions but it is not very accurate, I was hoping R would be more accurate and tell me what this line is, or at least a set of solutions. On Thu, May 6, 2021 at 5:28 AM Abbs Spurdle wrote: > Just realized five minutes after posting that I misinterpreted your > question, slightly. > However, after comparing the solution sets for *both* equations, I > can't see any obvious difference between the two. > If there is any difference, presumably that difference is extremely small. > > > On Thu, May 6, 2021 at 8:39 PM Abbs Spurdle wrote: > > > > Hi Ursula, > > > > If I'm not mistaken, there's an infinite number of solutions, which > > form a straight (or near straight) line. > > Refer to the following code, and attached plot. > > > > begin code--- > > library (barsurf) > > vF1 <- function (u, v) > > { n <- length (u) > > k <- numeric (n) > > for (i in seq_len (n) ) > > k [i] <- intfun1 (c (u [i], v [i]) ) > > k > > } > > plotf_cfield (vF1, c (0, 0.2), fb = (-2:2) / 10, > > main="(integral_1 - 1)", > > xlab="S[1]", ylab="S[2]", > > n=40, raster=TRUE, theme="heat", contour.labels=TRUE) > > end code > > > > I'm not familiar with the RootSolve package. > > Nor am I quite sure what you're trying to compute, given the apparent > > infinite set of solutions. > > > > So, for now at least, I'll leave comments on the root finding to someone > who is. > > > > > > Abby > > > > > > On Thu, May 6, 2021 at 8:46 AM Ursula Trigos-Raczkowski > > wrote: > > > > > > Hello, > > > I am trying to solve a system of integral equations using multiroot. I > have > > > tried asking on stack exchange and reddit without any luck. > > > Multiroot uses the library(RootSolve). > > > > > > I have two integral equations involving constants S[1] and S[2] (which > are > > > free.) I would like to find what *positive* values of S[1] and S[2] > make > > > the resulting > > > (Integrals-1) = 0. > > > (I know that the way I have the parameters set up the equations are > very > > > similar but I am interested in changing the parameters once I have the > code > > > working.) > > > My attempt at code: > > > > > > ```{r} > > > a11 <- 1 #alpha_{11} > > > a12 <- 1 #alpha_{12} > > > a21 <- 1 #alpha_{21} > > > a22 <- 1 #alpha_{22} > > > b1 <- 2 #beta1 > > > b2 <- 2 #beta2 > > > d1 <- 1 #delta1 > > > d2 <- 1 #delta2 > > > g <- 0.5 #gamma > > > > > > > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > > > > > #defining equation we would like to solve > > > intfun1<- function(S) {integrate(function(x) integrand1(x, > > > S),lower=0,upper=Inf)[[1]]-1} > > > intfun2<- function(S) {integrate(function(x) integrand2(x, > > > S),lower=0,upper=Inf)[[1]]-1} > > > > > > #putting both equations into one term > > > model <- function(S) c(F1 = intfun1,F2 = intfun2) > > > > > > #Solving for roots > > > (ss <-multiroot(f=model, start=c(0,0))) > > > ``` > > > > > > This gives me the error Error in stode(y, times, func, parms = parms, > ...) : > > > REAL() can only be applied to a 'numeric', not a 'list' > > > > > > However this simpler example works fine: > > > > > > ```{r} > > > #Defining the functions > > > model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) > > > > > > #Solving for the roots > > > (ss <- multiroot(f = model, start = c(0,0))) > > > ``` > > > > > > Giving me the required x_1= 4 and x_2 =1. > > > > > > I was given some code to perform a least squares analysis on the same > > > system but I neither understand the code, nor believe that it is doing > what > > > I am looking for as different initial values give wildly different S > values. > > > > > > ```{r} > > > a11 <- 1 #alpha_{11} > > > a12 <- 1 #alpha_{12} > > > a21 <- 1 #alpha_{21} > > > a22 <- 1 #alpha_{22} > > > b1 <- 2 #beta1 > > > b2 <- 2 #beta2 > > > d1 <- 1 #delta1 > > > d2 <- 1 #delta2 > > > g <- 0.5 #gamma > > > > > > > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > > > > > #defining equation we would like to solve > > > intfun1<- function(S) {integrate(function(x)integrand1(x, > > > S),lower=0,upper=Inf)[[1]]-1} > > > intfun2<- function(S) {integrate(function(x)integrand2(x, > > > S),lower=0,upper=Inf)[[1]]-1} > > > > > > #putting both equations into one term > > > model <- function(S) if(any(S<0))NA else intfun1(S)**2+ intfun2(S)**2 > > > > >
Re: [R] solving integral equations with undefined parameters using multiroot
Just realized five minutes after posting that I misinterpreted your question, slightly. However, after comparing the solution sets for *both* equations, I can't see any obvious difference between the two. If there is any difference, presumably that difference is extremely small. On Thu, May 6, 2021 at 8:39 PM Abbs Spurdle wrote: > > Hi Ursula, > > If I'm not mistaken, there's an infinite number of solutions, which > form a straight (or near straight) line. > Refer to the following code, and attached plot. > > begin code--- > library (barsurf) > vF1 <- function (u, v) > { n <- length (u) > k <- numeric (n) > for (i in seq_len (n) ) > k [i] <- intfun1 (c (u [i], v [i]) ) > k > } > plotf_cfield (vF1, c (0, 0.2), fb = (-2:2) / 10, > main="(integral_1 - 1)", > xlab="S[1]", ylab="S[2]", > n=40, raster=TRUE, theme="heat", contour.labels=TRUE) > end code > > I'm not familiar with the RootSolve package. > Nor am I quite sure what you're trying to compute, given the apparent > infinite set of solutions. > > So, for now at least, I'll leave comments on the root finding to someone who > is. > > > Abby > > > On Thu, May 6, 2021 at 8:46 AM Ursula Trigos-Raczkowski > wrote: > > > > Hello, > > I am trying to solve a system of integral equations using multiroot. I have > > tried asking on stack exchange and reddit without any luck. > > Multiroot uses the library(RootSolve). > > > > I have two integral equations involving constants S[1] and S[2] (which are > > free.) I would like to find what *positive* values of S[1] and S[2] make > > the resulting > > (Integrals-1) = 0. > > (I know that the way I have the parameters set up the equations are very > > similar but I am interested in changing the parameters once I have the code > > working.) > > My attempt at code: > > > > ```{r} > > a11 <- 1 #alpha_{11} > > a12 <- 1 #alpha_{12} > > a21 <- 1 #alpha_{21} > > a22 <- 1 #alpha_{22} > > b1 <- 2 #beta1 > > b2 <- 2 #beta2 > > d1 <- 1 #delta1 > > d2 <- 1 #delta2 > > g <- 0.5 #gamma > > > > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > > > #defining equation we would like to solve > > intfun1<- function(S) {integrate(function(x) integrand1(x, > > S),lower=0,upper=Inf)[[1]]-1} > > intfun2<- function(S) {integrate(function(x) integrand2(x, > > S),lower=0,upper=Inf)[[1]]-1} > > > > #putting both equations into one term > > model <- function(S) c(F1 = intfun1,F2 = intfun2) > > > > #Solving for roots > > (ss <-multiroot(f=model, start=c(0,0))) > > ``` > > > > This gives me the error Error in stode(y, times, func, parms = parms, ...) : > > REAL() can only be applied to a 'numeric', not a 'list' > > > > However this simpler example works fine: > > > > ```{r} > > #Defining the functions > > model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) > > > > #Solving for the roots > > (ss <- multiroot(f = model, start = c(0,0))) > > ``` > > > > Giving me the required x_1= 4 and x_2 =1. > > > > I was given some code to perform a least squares analysis on the same > > system but I neither understand the code, nor believe that it is doing what > > I am looking for as different initial values give wildly different S values. > > > > ```{r} > > a11 <- 1 #alpha_{11} > > a12 <- 1 #alpha_{12} > > a21 <- 1 #alpha_{21} > > a22 <- 1 #alpha_{22} > > b1 <- 2 #beta1 > > b2 <- 2 #beta2 > > d1 <- 1 #delta1 > > d2 <- 1 #delta2 > > g <- 0.5 #gamma > > > > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > > > #defining equation we would like to solve > > intfun1<- function(S) {integrate(function(x)integrand1(x, > > S),lower=0,upper=Inf)[[1]]-1} > > intfun2<- function(S) {integrate(function(x)integrand2(x, > > S),lower=0,upper=Inf)[[1]]-1} > > > > #putting both equations into one term > > model <- function(S) if(any(S<0))NA else intfun1(S)**2+ intfun2(S)**2 > > > > #Solving for roots > > optim(c(0,0), model) > > ``` > > > > I appreciate any tips/help as I have been struggling with this for some > > weeks now. > > thank you, > > -- > > Ursula > > Ph.D. student, University of Michigan > > Applied and Interdisciplinary Mathematics > > utri...@umich.edu > > > > [[alternative HTML version deleted]] > > > > __ > > R-help@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.
Re: [R] solving integral equations with undefined parameters using multiroot
Hi Ursula, If I'm not mistaken, there's an infinite number of solutions, which form a straight (or near straight) line. Refer to the following code, and attached plot. begin code--- library (barsurf) vF1 <- function (u, v) { n <- length (u) k <- numeric (n) for (i in seq_len (n) ) k [i] <- intfun1 (c (u [i], v [i]) ) k } plotf_cfield (vF1, c (0, 0.2), fb = (-2:2) / 10, main="(integral_1 - 1)", xlab="S[1]", ylab="S[2]", n=40, raster=TRUE, theme="heat", contour.labels=TRUE) end code I'm not familiar with the RootSolve package. Nor am I quite sure what you're trying to compute, given the apparent infinite set of solutions. So, for now at least, I'll leave comments on the root finding to someone who is. Abby On Thu, May 6, 2021 at 8:46 AM Ursula Trigos-Raczkowski wrote: > > Hello, > I am trying to solve a system of integral equations using multiroot. I have > tried asking on stack exchange and reddit without any luck. > Multiroot uses the library(RootSolve). > > I have two integral equations involving constants S[1] and S[2] (which are > free.) I would like to find what *positive* values of S[1] and S[2] make > the resulting > (Integrals-1) = 0. > (I know that the way I have the parameters set up the equations are very > similar but I am interested in changing the parameters once I have the code > working.) > My attempt at code: > > ```{r} > a11 <- 1 #alpha_{11} > a12 <- 1 #alpha_{12} > a21 <- 1 #alpha_{21} > a22 <- 1 #alpha_{22} > b1 <- 2 #beta1 > b2 <- 2 #beta2 > d1 <- 1 #delta1 > d2 <- 1 #delta2 > g <- 0.5 #gamma > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > #defining equation we would like to solve > intfun1<- function(S) {integrate(function(x) integrand1(x, > S),lower=0,upper=Inf)[[1]]-1} > intfun2<- function(S) {integrate(function(x) integrand2(x, > S),lower=0,upper=Inf)[[1]]-1} > > #putting both equations into one term > model <- function(S) c(F1 = intfun1,F2 = intfun2) > > #Solving for roots > (ss <-multiroot(f=model, start=c(0,0))) > ``` > > This gives me the error Error in stode(y, times, func, parms = parms, ...) : > REAL() can only be applied to a 'numeric', not a 'list' > > However this simpler example works fine: > > ```{r} > #Defining the functions > model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) > > #Solving for the roots > (ss <- multiroot(f = model, start = c(0,0))) > ``` > > Giving me the required x_1= 4 and x_2 =1. > > I was given some code to perform a least squares analysis on the same > system but I neither understand the code, nor believe that it is doing what > I am looking for as different initial values give wildly different S values. > > ```{r} > a11 <- 1 #alpha_{11} > a12 <- 1 #alpha_{12} > a21 <- 1 #alpha_{21} > a22 <- 1 #alpha_{22} > b1 <- 2 #beta1 > b2 <- 2 #beta2 > d1 <- 1 #delta1 > d2 <- 1 #delta2 > g <- 0.5 #gamma > > > integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* > x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} > integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* > x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} > > #defining equation we would like to solve > intfun1<- function(S) {integrate(function(x)integrand1(x, > S),lower=0,upper=Inf)[[1]]-1} > intfun2<- function(S) {integrate(function(x)integrand2(x, > S),lower=0,upper=Inf)[[1]]-1} > > #putting both equations into one term > model <- function(S) if(any(S<0))NA else intfun1(S)**2+ intfun2(S)**2 > > #Solving for roots > optim(c(0,0), model) > ``` > > I appreciate any tips/help as I have been struggling with this for some > weeks now. > thank you, > -- > Ursula > Ph.D. student, University of Michigan > Applied and Interdisciplinary Mathematics > utri...@umich.edu > > [[alternative HTML version deleted]] > > __ > R-help@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. __ R-help@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.
[R] solving integral equations with undefined parameters using multiroot
Hello, I am trying to solve a system of integral equations using multiroot. I have tried asking on stack exchange and reddit without any luck. Multiroot uses the library(RootSolve). I have two integral equations involving constants S[1] and S[2] (which are free.) I would like to find what *positive* values of S[1] and S[2] make the resulting (Integrals-1) = 0. (I know that the way I have the parameters set up the equations are very similar but I am interested in changing the parameters once I have the code working.) My attempt at code: ```{r} a11 <- 1 #alpha_{11} a12 <- 1 #alpha_{12} a21 <- 1 #alpha_{21} a22 <- 1 #alpha_{22} b1 <- 2 #beta1 b2 <- 2 #beta2 d1 <- 1 #delta1 d2 <- 1 #delta2 g <- 0.5 #gamma integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} #defining equation we would like to solve intfun1<- function(S) {integrate(function(x) integrand1(x, S),lower=0,upper=Inf)[[1]]-1} intfun2<- function(S) {integrate(function(x) integrand2(x, S),lower=0,upper=Inf)[[1]]-1} #putting both equations into one term model <- function(S) c(F1 = intfun1,F2 = intfun2) #Solving for roots (ss <-multiroot(f=model, start=c(0,0))) ``` This gives me the error Error in stode(y, times, func, parms = parms, ...) : REAL() can only be applied to a 'numeric', not a 'list' However this simpler example works fine: ```{r} #Defining the functions model <- function(x) c(F1 = x[1]+ 4*x[2] -8,F2 = x[1]-4*x[2]) #Solving for the roots (ss <- multiroot(f = model, start = c(0,0))) ``` Giving me the required x_1= 4 and x_2 =1. I was given some code to perform a least squares analysis on the same system but I neither understand the code, nor believe that it is doing what I am looking for as different initial values give wildly different S values. ```{r} a11 <- 1 #alpha_{11} a12 <- 1 #alpha_{12} a21 <- 1 #alpha_{21} a22 <- 1 #alpha_{22} b1 <- 2 #beta1 b2 <- 2 #beta2 d1 <- 1 #delta1 d2 <- 1 #delta2 g <- 0.5 #gamma integrand1 <- function(x,S) {b1*g/d1*exp(-g*x)*(1-exp(-d1* x))*exp(-a11*b1*S[1]/d1*(1-exp(-d1*x))-a12*b2*S[2]/d2*(1-exp(-d2*x)))} integrand2 <- function(x,S) {b2*g/d2*exp(-g*x)*(1-exp(-d2* x))*exp(-a22*b2*S[2]/d2*(1-exp(-d2*x))-a21*b1*S[1]/d1*(1-exp(-d1*x)))} #defining equation we would like to solve intfun1<- function(S) {integrate(function(x)integrand1(x, S),lower=0,upper=Inf)[[1]]-1} intfun2<- function(S) {integrate(function(x)integrand2(x, S),lower=0,upper=Inf)[[1]]-1} #putting both equations into one term model <- function(S) if(any(S<0))NA else intfun1(S)**2+ intfun2(S)**2 #Solving for roots optim(c(0,0), model) ``` I appreciate any tips/help as I have been struggling with this for some weeks now. thank you, -- Ursula Ph.D. student, University of Michigan Applied and Interdisciplinary Mathematics utri...@umich.edu [[alternative HTML version deleted]] __ R-help@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.