1

If I have a function 

estimator <- function(A,B) {
 A*(B+23)
}

How can I reverse this function to find the value of A for B as a sequence between 0 and 120 (B=1,2,3,4,...,120) that would give a fixed result, say C = 20?

I would use it to map the values for which satisfy the equation A*(B+23)= C = 20 with B being a list b.list between 0 and 120, for c.list, of different C?

b.list <- seq(0,120,by=1)
c.list <- tibble(seq(10,32,by=2))

In the end, I would like to plot the lines of curves of the function for different C using purrr or similar.

I.e.: Given that the height of a tree in metres at age 100 will follow the function, C = A*(B+23), solve for A that will give the result C=10 when B, Age is a list of years between 0 and 120?

Here's a link showing what I'm trying to make! enter image description here

Here's another one

enter image description here Many thanks!

Silviculturalist
  • 199
  • 11

1 Answers1

0

For the inverse it is a quick inversion :

A = C/(B+23)

One answer could be :

B <- seq(0, 120)
C <- seq(10, 32, 2)
A <- matrix(0,
            nrow = length(B),
            ncol = length(C))

for(i in 1:ncol(M)){
    A[,i] <- C[i] / (B + 23)
}

matplot(B, A, type ="l", col = "black")

In case of a more complex function indeed you need an automatic solving problem. One way is to see it like an optimisation problem where you want to minimise the distance from C :

B <- seq(1, 120)
C <- seq(10, 32, 2)
A <- matrix(0,
            nrow = length(B),
            ncol = length(C))
fct <- function(A, B, C){
    paramasi <- 25
    parambeta<- 7395.6
    paramb2 <- -1.7829
    refB <- 100
    d <- parambeta*(paramasi^paramb2)
    r <- (((A-d)^2)+(4*parambeta*A*(B^paramb2)))^0.5
    si_est <- (A+d+r)/ (2+(4*parambeta*(refB^paramb2)) / (A-d+r))
    return(sum(si_est - C)^2)}

for(c in 1:length(C)){
    for(b in 1:length(B)){
        # fixe parameters + optimisation
        res <- optim(par = 1, fn = fct, B = B[b], C = C[c])
        A[b, c] <- res$par
    }
}

matplot(B, A, type = "l", col = "black")

You need to be careful because in your case I think that you could find an analytical formula for the inverse which would be better.

Good luck !

Rémi Coulaud
  • 1,684
  • 1
  • 8
  • 19
  • Well yes, this inverts the formula, but not the function! A is not really a parameter. The real formula would be the inverse of a function for one species, approximation of the site index (the height of a tree at a reference age, e.g. 100) ?```function(topheight, age){ paramasi <- 25 parambeta <- 7395.6 paramb2 <- -1.7829 refAge <- 100 d <- parambeta*(paramasi^paramb2) r <- (((topheight-d)^2)+(4*parambeta*topheight*(age^paramb2)))^0.5 si_est<- (topheight+d+r)/ (2+(4*parambeta*(refAge^paramb2)) / (topheight-d+r)) return(si_est) }``` – Silviculturalist May 06 '20 at 12:00
  • Be careful the result of your function need to be a vector whereas you only give a result under the form of a constant : 20 – Rémi Coulaud May 06 '20 at 12:33