R:尝试编写等效函数时出错 n 选择 k
R: error when trying to write an equivalent function n choose k
我正在学习 class R 编程入门。
我们被要求编写一个与 n choose k:
相同的函数
choose(n, k)
我们被要求检查函数是否工作 运行 n = 200,k = 50。
我写了下面的代码:
select_k <- function(n, k){
sr <- c(log10(log10((factorial(n-1)/factorial(k-1)*factorial(n-k-2)))*(n/k)))
return(sr)
}
因为select_k应该是“n选k”。
我的函数适用于以下值:100 选择 25,但不适用于更大的值,例如 n = 200、k = = 50。
select_k( n = 200, k = 50)
[1] NaN
Warning message:
In factorial(n) : value out of range in 'gammafn'
我不知道还能做些什么来解决这个问题。
这不适用于较大的 n
,因为 factorial(n)
太大了:
> factorial(199)
[1] Inf
Warning message:
In factorial(199) : value out of range in 'gammafn'
这应该是 return 200,但计算机只看到您试图将 Inf
除以 Inf
:
> factorial(200)/factorial(199)
[1] NaN
Warning messages:
1: In factorial(200) : value out of range in 'gammafn'
2: In factorial(199) : value out of range in 'gammafn'
显然 "n choose k" 中的很多乘法抵消了,因此您需要避免使用常规阶乘,而只乘以不抵消的数字(?prod
可能会有用为你)。或者(可能更好)使用日志版本 lfactorial
来避免 运行 变成您的计算机无法存储的数字。
编辑:添加了来自@MrFlick 评论的lfactorial
推荐
看看这个{
a <- function(n, k) {
exp(lgamma(n+1) - lgamma(n - k + 1) - lgamma(k + 1) )
}
我正在学习 class R 编程入门。
我们被要求编写一个与 n choose k:
相同的函数choose(n, k)
我们被要求检查函数是否工作 运行 n = 200,k = 50。
我写了下面的代码:
select_k <- function(n, k){
sr <- c(log10(log10((factorial(n-1)/factorial(k-1)*factorial(n-k-2)))*(n/k)))
return(sr)
}
因为select_k应该是“n选k”。
我的函数适用于以下值:100 选择 25,但不适用于更大的值,例如 n = 200、k = = 50。
select_k( n = 200, k = 50)
[1] NaN
Warning message:
In factorial(n) : value out of range in 'gammafn'
我不知道还能做些什么来解决这个问题。
这不适用于较大的 n
,因为 factorial(n)
太大了:
> factorial(199)
[1] Inf
Warning message:
In factorial(199) : value out of range in 'gammafn'
这应该是 return 200,但计算机只看到您试图将 Inf
除以 Inf
:
> factorial(200)/factorial(199)
[1] NaN
Warning messages:
1: In factorial(200) : value out of range in 'gammafn'
2: In factorial(199) : value out of range in 'gammafn'
显然 "n choose k" 中的很多乘法抵消了,因此您需要避免使用常规阶乘,而只乘以不抵消的数字(?prod
可能会有用为你)。或者(可能更好)使用日志版本 lfactorial
来避免 运行 变成您的计算机无法存储的数字。
编辑:添加了来自@MrFlick 评论的lfactorial
推荐
看看这个{
a <- function(n, k) {
exp(lgamma(n+1) - lgamma(n - k + 1) - lgamma(k + 1) )
}