Implementing an algebraic equation in Base R

46 views Asked by At

I'm new to algebraic equations in R and am trying to implement the following equation in R (below).

However, as precise as I have tried to be, I expect the answer to be 0.050865 but instead I get 11.43653.

I wonder what I am missing as I really can't find any mistakes?

delta_T=0.1522
n=18
NT=324
NC=162
p=0.264
N=NT+NC

V_T2 =  ((NT+NC)/(NT*NC))*(1 + ((n - 1)*p)) +
    ( delta_T^2 * ( 
      (((N - 2)*(1-p)^2 ) + (n*(N - 2*n)*p^2) +
         (2* (N - 2*n) * p * (1 - p)) ) /
         (2* (N-2)) * ( (N-2) - (2* (n-1)*p) ) 
    )  )

V_T2

# > [1] 11.43653 BUT I expect the answer to be 0.050865

enter image description here

1

There are 1 answers

0
Mark On BEST ANSWER

The reason you're getting the wrong answer is because, in the second term, you think you're doing

numerator 
/
denominator

but what you're actually doing is

numerator 
/
(2 * (N - 2)) 
* 
((N - 2) - (2 * (n - 1) * p) )

So, in other words, you were missing parentheses around the denominator.

But, in general, it's helpful to break your code down into intermediate variables, when dealing with complicated equations, so it's easier to debug, and you can do sanity checks on them:

equation <- function(NT, NC, n, p, delta_T) {
    N <- NT + NC
    first_term <- N / (NT * NC) * (1 + ((n - 1) * p))
    numerator <- (N - 2) * (1 - p)^2 + n * (N - 2 * n) * p^2 +
         2 * (N - 2 * n) * p * (1 - p)
    denominator <- 2 * (N-2) * (N - 2 - 2 * (n - 1) * p) 
    second_term <- delta_T^2 * (numerator / denominator)
    first_term + second_term
}