2

我将从我已经做过的开始。我正在寻找一种通过更改每行的b和s参数来求解方程f的方法， Q和n是常数。我知道 apply() 适用于此类问题，但这似乎对我不起作用。我想找到的变量没有给出唯一的解决方案。

Q = 0.203
n = 0.014
f <- function(y) (Q - (1/n)*(y*b)*((y*b)/(2*y+b))^(2/3)*sqrt(s))

使用这些参数，假设b = 0.5并s = 0.01使用 uniroot() 我得到以下信息。这就是我想要的结果。

uniroot(f, lower = 0.000001, upper = 1000000)$root

[1] 0.2328931

（那些较低和较高的值似乎对我很有效）

现在我需要为大型数据集解决此功能。

set.seed(123)
tibble::tibble(b = runif(n = 1000, min = 0.1, max = 1.5),
               s = runif(n = 1000, min = 0.001, max = 5)) %>% 
  dplyr::mutate(yn = uniroot(f, lower = 0.000001, upper = 1000000)$root) %>% 
  head(5)

这是我想要的输出。

      b     s    yn
1 0.503 1.37  0.0434
2 1.20  2.97  0.0194
3 0.673 0.802 0.0421
4 1.34  4.27  0.0163
5 1.42  4.24  0.0157
4

2 回答 2

2

考虑更改函数以从数据中获取“b”、“s”列并使用rowwise

f <- function(y, dat) with(dat, (Q - (1/n)*(y*b)*((y*b)/(2*y+ b))^(2/3)*sqrt(s)))
df1 %>%
    rowwise %>%
    dplyr::mutate(yn = uniroot(f, lower = 0.000001, upper = 1000000, 
            dat = cur_data())$root) %>% 
   ungroup %>%
   head(5)

-输出

# A tibble: 5 x 3
#      b     s     yn
#  <dbl> <dbl>  <dbl>
#1 0.503 1.37  0.0434
#2 1.20  2.97  0.0194
#3 0.673 0.802 0.0421
#4 1.34  4.27  0.0163
#5 1.42  4.24  0.0157

或者另一个选项pmap来自purrr

library(purrr)
df1 %>%
  mutate(yn = pmap_dbl(select(., b, s), ~
       uniroot(f,  lower = 0.000001, upper = 1000000,
        dat = tibble(b = ..1, s = ..2))$root))

-输出

# A tibble: 10,000 x 3
#       b     s     yn
#   <dbl> <dbl>  <dbl>
# 1 0.503 1.37  0.0434
# 2 1.20  2.97  0.0194
# 3 0.673 0.802 0.0421
# 4 1.34  4.27  0.0163
# 5 1.42  4.24  0.0157
# ...

数据

set.seed(123)
df1 <- tibble::tibble(b = runif(n = 1000, min = 0.1, max = 1.5),
               s = runif(n = 1000, min = 0.001, max = 5))
于 2020-12-08T22:22:48.767 回答
2

实际上，您已经非常接近目标了。这是使用Vectorize+的基本 R 选项do.call，它可以帮助您

f <- function(b, s) {
  fn <- function(y) (Q - (1 / n) * (y * b) * ((y * b) / (2 * y + b))^(2 / 3) * sqrt(s))
  uniroot(fn, lower = 0.000001, upper = 1000000)$root
}

df$yn <- do.call(Vectorize(f), df)

这样

> df
# A tibble: 1,000 x 3
       b     s     yn
   <dbl> <dbl>  <dbl>
 1 0.503 1.37  0.0435
 2 1.20  2.97  0.0194
 3 0.673 0.802 0.0422
 4 1.34  4.27  0.0163
 5 1.42  4.24  0.0157
 6 0.164 2.39  0.0912
 7 0.839 3.87  0.0224
 8 1.35  1.48  0.0223
 9 0.872 0.329 0.0468
10 0.739 2.20  0.0289
# ... with 990 more rows

数据

set.seed(123)
df <- tibble::tibble(b = runif(n = 1000, min = 0.1, max = 1.5),
               s = runif(n = 1000, min = 0.001, max = 5))
于 2020-12-08T23:36:56.447 回答