[R] R emulation of FindRoot in Mathematica

Thu Jan 19 18:10:15 CET 2023

Dear Troels,

Are you aware of the caracas package for computer algebra in R? The package provides an interface to the SymPy package in python via the reticulate package.

I have no idea if the package can be helpful in your connection, but please report succeses and failures back.

Best regards
Søren

-----Original Message-----
From: Troels Ring <tring using gvdnet.dk<mailto:Troels%20Ring%20%3ctring using gvdnet.dk%3e>>
To: Jeff Newmiller <jdnewmil using dcn.davis.ca.us<mailto:Jeff%20Newmiller%20%3cjdnewmil using dcn.davis.ca.us%3e>>, r-help using r-project.org<mailto:r-help using r-project.org>, "Ebert,Timothy Aaron" <tebert using ufl.edu<mailto:%22Ebert,Timothy%20Aaron%22%20%3ctebert using ufl.edu%3e>>, Valentin Petzel <valentin using petzel.at<mailto:Valentin%20Petzel%20%3cvalentin using petzel.at%3e>>
Subject: Re: [R] R emulation of FindRoot in Mathematica
Date: Thu, 19 Jan 2023 16:18:57 +0100

Hi Jeff - that is definitely not fair, this is a highly respected

scientist but old me are probably not clever enough to explain problem

which I think is not too difficult. Basically, it is a series of

consecutive statements of H and Mg and K combining to ATP and ADP,

creatin and creatinP under assumption og known total inorganic phosphate

(pi), total ATP total ADP, creatin and creatinphosphate. It was written

in 1997 when R was less mature - and if we know how to do it in R so not

to depend on black box operations it would be fun. For the moment, I'm

not even sure I understand why there is now focus on the problem being

nonlinear, since I want to evaluate at given H or pH.  hatp <- 10^6.494

* H * ATP ?

BW Troels

Den 19-01-2023 kl. 15:54 skrev Jeff Newmiller:

But it is simultaneously an example of why some researchers like black box solvers... a system of dozens of nonlinear equations can potentially have many or even infinite solutions. If the researcher is weak in math, they may have no idea which solutions are possible and having a tool like FindRoot confidently return a solution lets them focus on other things. Sort of like ChatGPT.

TL;DR the author may have no idea about how to resolve this without relying on the opaque FindRoot.

On January 19, 2023 6:28:53 AM PST, "Ebert,Timothy Aaron" <

<mailto:tebert using ufl.edu>

tebert using ufl.edu

> wrote:

This is a poster child for why we like open source software. "I dump numbers into a black box and get numbers out but I cannot verify how the numbers out were calculated so they must be correct" approach to analysis does not really work for me.

Tim

-----Original Message-----

From: R-help <

<mailto:r-help-bounces using r-project.org>

r-help-bounces using r-project.org

> On Behalf Of Troels Ring

Sent: Thursday, January 19, 2023 9:18 AM

To: Valentin Petzel <

<mailto:valentin using petzel.at>

valentin using petzel.at

>; r-help mailing list <

<mailto:r-help using r-project.org>

r-help using r-project.org

>

Subject: Re: [R] R emulation of FindRoot in Mathematica

[External Email]

Thanks,   Valentin for the suggestion. I'm not sure I can go that way. I

include below the statements from the paper containing the knowledge on the basis of which I would like to know at specified [H] the concentration of each of the many metabolites given the constraints. I have tried to contact the author to get the full code but it seems difficult.

BW Troels

hatp <- 10^6.494*H*atp

hhatp <- 10^3.944*H*hatp

hhhatp <- 10^1.9*H*hhatp

hhhhatp  <- 10*H*hhhatp

mgatp <- 10^4.363*atp*mg

mghatp <- 10^2.299*hatp*mg

mg2atp <- 10^1-7*mg*mgatp

katp <- 10^0.959*atp*k

hadp <- 10^6.349*adp*H

hhadp <- 10^3.819*hadp*H

hhhadp <- 10*H*hhadp

mgadp <- 10^3.294*mg*adp

mghadp <- 10^1.61*mg*hadp

mg2adp <- 10*mg*mgadp

kadp <- 10^0.82*k*adp

hpi <- 10^11.616*H*pi

hhpi <- 10^6.7*h*hpi

hhhpi <- 10^1.962*h*hhpi

mgpi <- 10^3.4*mg*pi

mghpi <- 10^1.946*mg*hpi

mghhpi <- 10^1.19*mg*hhpi

kpi <- 10^0.6*k*pi

khpi <- 10^1.218*k*hpi

khhpi <- 10^-0.2*k*hhpi

hpcr <- 10^14.3*h*pcr

hhpcr <- 10^4.5*h*hpcr

hhhpcr <- 10^2.7*h*hhpcr

hhhhpcr <- 100*h*hhhpcr

mghpcr <- 10^1.6*mg*hpcr

kpcr <- 10^0.74*k*pcr

khpcr <- 10^0.31*k*hpcr

khhpcr <- 10^-0.13*k*hhpcr

hcr <- 10^14.3*h*cr

hhcr <- 10^2.512*h*hcr

hlactate <- 10^3.66*h*lactate

mglactate <- 10^0.93*mg*lactate

tatp <- atp + hatp + hhatp + hhhatp + mgatp + mghatp + mg2atp + katp

tadp <- adp + hadp + hhadp + hhhadp + mghadp + mgadp + mg2adp + kadp

tpi <- pi + hpi + hhpi + hhhpi + mgpi + mghpi + mghhpi + kpi + khpi + khhpi

tpcr <- pcr + hpcr + hhpcr + hhhpcr + hhhhpcr + mghpcr + kpcr + khpcr + khhpcr

tcr <- cr + hcr + hhcr

tmg <- mg + mgatp + mghatp + mg2atp + mgadp + mghadp + mg2adp + mgpi + kghpi + mghhpi +

   mghpcr + mglactate

tk <- k + katp + kadp + kpi + khpi + khhpi + kpcr + khpcr + khhpcr

tlactate <- lactate + hlactate + mglactate

# conditions

tatp <- 0.008

tpcr <- 0.042

tcr <- 0.004

tadp <- 0.00001

tpi <- 0.003

tlactate <- 0.005

# free K and Mg constrained to be fixed

#

mg <- 0.0006

k <- 0.12

Den 19-01-2023 kl. 12:11 skrev Valentin Petzel:

Hello Troels,

As fair as I understand you attempt to numerically solve a system of

non linear equations in multiple variables in R. R does not provide

this functionality natively, but have you tried multiroot from the

rootSolve package:

<https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fcran>

https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fcran

.r-project.org%2Fweb%2Fpackages%2FrootSolve%2FrootSolve.pdf&data=05%7C

01%7Ctebert%40ufl.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a3

14d76ace60a62331e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZ

sb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3

D%7C3000%7C%7C%7C&sdata=D9A3fwJ5x7GbEV4A01wncLUil7szTdSPul5vd0lsSBw%3D

&reserved=0

multiroot is called like

multiroot(f, start, ...)

where f is a function of one argument which is a vector of n values

(representing the n variables) and returning a vector of d values

(symbolising the d equations) and start is a vector of length n.

E.g. if we want so solve

x^2 + y^2 + z^2 = 1

x^3-y^3 = 0

x - z = 0

(which is of course equivalent to x = y = z, x^2 + y^2 + z^2 = 1, so x

= y = z = ±sqrt(1/3) ~ 0.577)

we'd enter

f <- function(x) c(x[1]**2 + x[2]**2 + x[3]**2 - 1, x[1]**3 - x[2]**3,

x[1] - x[3])

multiroot(f, c(0,0,0))

which yields

$root

[1] 0.5773502 0.5773505 0.5773502

$f.root

[1] 1.412261e-07 -2.197939e-07  0.000000e+00

$iter

[1] 31

$estim.precis

[1] 1.2034e-07

Best regards,

Valentin

Am Donnerstag, 19. Jänner 2023, 10:41:22 CET schrieb Troels Ring:

Hi friends - I hope this is not a misplaced question. From the

literature (Kushmerick AJP 1997;272:C1739-C1747) I have a series of

Mathematica equations which are solved together to yield over

different

pH values the concentrations of metabolites in skeletal muscle using

the

Mathematica function FindRoot((E1,E2...),(V2,V2..)] where E is a

list of

equations and V list of variables.  Most of the equations are

individual

binding reactions of the form 10^6.494*atp*h == hatp and next

10^9.944*hatp*h ==hhatp describing binding of singe protons or Mg or

K

to ATP or creatin for example, but we also have constraints giving

total

concentrations of say ATP i.e. ATP + ATPH, ATPH2..ATP.Mg

I have, without success, tried to find ways to do this in R - I have

36

equations on 36 variables and 8 equations on total concentrations.

As

far as I can see from the definition of FindRoot in Wolfram, Newton

search or secant search is employed.

I'm on Windows R 4.2.2

Best wishes

Troels Ring, MD

Aalborg, Denmark

______________________________________________

<mailto:R-help using r-project.org>

R-help using r-project.org

 mailing list -- To UNSUBSCRIBE and more, see

<https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fst>

https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fst

at.ethz.ch%2Fmailman%2Flistinfo%2Fr-help&data=05%7C01%7Ctebert%40ufl

.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a6233

1e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoi

MC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C

%7C%7C&sdata=7DTBQItdQpAqK%2FCS1%2BqQvYdlvjyJMjzTOXhoS6AY%2FJQ%3D&re

served=0

PLEASE do read the posting guide

<https://nam10.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.r>

https://nam10.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.r

-project.org%2Fposting-guide.html&data=05%7C01%7Ctebert%40ufl.edu%7C7c

b98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a62331e1b84%7C0%

7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiL

CJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=m0

0vZP75nk9icL6H8Gc0dH1bHhkRCS9I5N27uORQmQ0%3D&reserved=0

and provide commented, minimal, self-contained, reproducible code.

        [[alternative HTML version deleted]]

______________________________________________

<mailto:R-help using r-project.org>

R-help using r-project.org

 mailing list -- To UNSUBSCRIBE and more, see

<https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-help&data=05%7C01%7Ctebert%40ufl.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a62331e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=7DTBQItdQpAqK%2FCS1%2BqQvYdlvjyJMjzTOXhoS6AY%2FJQ%3D&reserved=0>

https://nam10.safelinks.protection.outlook.com/?url=https%3A%2F%2Fstat.ethz.ch%2Fmailman%2Flistinfo%2Fr-help&data=05%7C01%7Ctebert%40ufl.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a62331e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=7DTBQItdQpAqK%2FCS1%2BqQvYdlvjyJMjzTOXhoS6AY%2FJQ%3D&reserved=0

PLEASE do read the posting guide

<https://nam10.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.r-project.org%2Fposting-guide.html&data=05%7C01%7Ctebert%40ufl.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a62331e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=m00vZP75nk9icL6H8Gc0dH1bHhkRCS9I5N27uORQmQ0%3D&reserved=0>

https://nam10.safelinks.protection.outlook.com/?url=http%3A%2F%2Fwww.r-project.org%2Fposting-guide.html&data=05%7C01%7Ctebert%40ufl.edu%7C7cb98cd926b34284cd5f08dafa28026c%7C0d4da0f84a314d76ace60a62331e1b84%7C0%7C0%7C638097347110882622%7CUnknown%7CTWFpbGZsb3d8eyJWIjoiMC4wLjAwMDAiLCJQIjoiV2luMzIiLCJBTiI6Ik1haWwiLCJXVCI6Mn0%3D%7C3000%7C%7C%7C&sdata=m00vZP75nk9icL6H8Gc0dH1bHhkRCS9I5N27uORQmQ0%3D&reserved=0

and provide commented, minimal, self-contained, reproducible code.

______________________________________________

<mailto:R-help using r-project.org>

R-help using r-project.org

 mailing list -- To UNSUBSCRIBE and more, see

<https://stat.ethz.ch/mailman/listinfo/r-help>

https://stat.ethz.ch/mailman/listinfo/r-help

PLEASE do read the posting guide

<http://www.R-project.org/posting-guide.html>

http://www.R-project.org/posting-guide.html

and provide commented, minimal, self-contained, reproducible code.

______________________________________________

<mailto:R-help using r-project.org>

R-help using r-project.org

 mailing list -- To UNSUBSCRIBE and more, see

<https://stat.ethz.ch/mailman/listinfo/r-help>

https://stat.ethz.ch/mailman/listinfo/r-help

PLEASE do read the posting guide

<http://www.R-project.org/posting-guide.html>

http://www.R-project.org/posting-guide.html

and provide commented, minimal, self-contained, reproducible code.

--

Søren Højsgaard
Head of Department of Mathematical Sciences
Aalborg University, Denmark

	[[alternative HTML version deleted]]