Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to store output of a loop when resampling of parameter values for simultaneous equations in R?

I have a code that evaluates parameter values for 4 simultaneous equations. I'm particularly interested to store all the parameter combinations when a + b (stored in results$ab) is greater than 3000. If it is greater than 3000, then it's coded as "Yes". I want to write a for-loop that will loop through the code to check for if a + b > 3000 and store the corresponding values. Then, I want the program to loop 1000 times, and store the parameter values for corresponding "Yes". I am trying to store the output, but it's not yielding me any results.

x <- seq(from = 0.0001, to = 1000, by = 0.1)
t <- seq(from = 0.0001, to = 1000, by = 0.1)
v <- seq(from = 0.0001, to = 1000, by = 0.1)
w <- seq(from = 0.0001, to = 1000, by = 0.1)
n <- seq(from = 0.0001, to = 1000, by = 0.1)
f <- seq(from = 0.0001, to = 1000, by = 0.1)


  values <- list(x = x, t = t, v = v, w = w, n = n, f = f)
  for(i in 1:1000){
    eqs <- list(
      a = expression(x * t - 2 * x),
      b = expression(v - x^2), 
      c = expression(x - w*t - t*t), 
      d = expression((n - f)/t)
    )

    for(i in 1:1000){
    samples <- 10000
    values.sampled <- lapply(values, sample, samples)
    results[i] <- sapply(eqs, eval, envir = values.sampled)
    results[i] <- data.frame(results)
    results$ab[i] <- results$a[i] + results$b[i]
    results$binary[i] <- ifelse(results$ab[i] > 3000, "Yes","No")
    output[i] <- results[results$binary=="Yes",]

  }

what <- as.list(output)
like image 717
Biotechgeek Avatar asked Dec 07 '25 08:12

Biotechgeek


1 Answers

a+b is equal to (x * t - 2 * x) + (v - x^2), which is just a quadratic in x, so you can solve a+b>3000 analytically for x, v and t.

The inequality is x^2 + (2-t)x + (3000-v) < 0.

Substitute T = 2-t and V = 3000-v, then x^2 + Tx + V < 0.

For this to have any values less than zero, it must have two real roots, which implies that T^2 - 4V > 0 - i.e V < (T^2)/4. (https://en.wikipedia.org/wiki/Quadratic_formula)

Given T and V satisfying this inequality, the values of x for which a+b>3000 are those between the roots of the quadratic, i.e. |2x+T| < sqrt(T^2 - 4V).

So if you want a selection of values that satisfy the condition, it should be straightforward to loop through a range of values of t, pick a value of v that satisfies V < (T^2)/4, and then pick an x that falls in the appropriate range.

Here is one way...

t <- 1:1000
T <- 2 - t
V <- sapply((T ^ 2) / 4, function(z) runif(min = 0, max = z, n = 1)) #assumes V>0 (???)
v <- 3000 - V
x <- (sapply(sqrt(T ^ 2 - 4 * V), function(z) runif(min = -z, max = z, n = 1)) - T) / 2
ab <- (x * t - 2 * x) + (v - x ^ 2) #all >3000 (except for t=2, where ab=3000 exactly)
like image 119
Andrew Gustar Avatar answered Dec 08 '25 20:12

Andrew Gustar