0

看了这个问题后：R中的数值比较困难

我仍然卡住了，因为我依赖于在某些我无法编​​辑的函数中进行的等式比较（或者我可以吗？）

我在本地环境中测试三个数字的总和是否等于一个sum(p)==1 --> TRUE（功能到下一个 - 这可能吗？

更多细节：我正在尝试“优化”输入 CART 模型的先验概率，使用优化器 ( dfoptim package, nmkb) 选择先验组合，将它们发送到 rpart 包进行模型拟合，然后验证（rps函数）包得分 - 但在rpart包的某个地方，我的先验概率抛出了一个错误，因为rpart认为它们的总和不等于 1。

这是一个可重现的示例：

require('rpart')
require('verification')
require('dfoptim')
data(iris)
set.seed(1)
tmp1 <- paste0(names(iris),collapse="+")
tmp2 <- gsub("\\+Species","",tmp1)
fmlatext <- paste0("Species~",tmp2)
tree <- rpart(as.formula(fmlatext),data=iris,method="class")
objfun <- function(priors,fmlatext,data){
  p <- priors/sum(priors) # turn arbitrary threesome into numbers that sum to 1
  p[1] <- 1-(sum(p)-p[1]) # ensure that numbers sum to 1
  print(c(p,sum(p)),digits=16)
  tree <- rpart(as.formula(fmlatext),data=data,parms=list(prior=p),
                method="class") 
  rpst <- rps(data$Species,predict(tree,data=data))
  return(rpst$rpss)
}
nlev <- nlevels(iris$Species)
guess <- seq(nlev)*10
lb <- rep(1,nlev)
ub <- rep(100,nlev)
bestpriors <- nmkb(par=guess,fn=objfun,lower=lb,upper=ub,
                   control=list(maximize=TRUE),fmlatext=fmlatext,data=iris)

运行这段代码给了我这个输出：

[1] 0.1666666666666667 0.3333333333333333 0.5000000000000000 1.0000000000000000
[1] 0.4353687449261023 0.2354416940871099 0.3291895609867877 1.0000000000000000
[1] 0.1224920651311070 0.5548713793562775 0.3226365555126156 1.0000000000000000
[1] 0.1268712138061573 0.2390044736120877 0.6341243125817551 1.0000000000000000
[1] 0.35141687748184969 0.57028058689316308 0.07830253562498726 1.00000000000000000
[1] 0.2997590406445614 0.5077659444797995 0.1924750148756391 1.0000000000000000
[1] 0.3598141573675122 0.4350423262345758 0.2051435163979119 0.9999999999999999
Error in get(paste("rpart", method, sep = "."), envir = environment())(Y,  : 
  Priors must sum to 1

在我的真实代码中，这种情况不一致地发生，具体取决于数据和猜测值，但它确实发生了，而且真的很痛苦。

我怎样才能克服这个错误？干杯，R

4

1 回答 1

0

潜在答案；不确定它是否健壮，但我创建了这个函数，它适用于一些以前对我来说失败的测试用例：

makeSumToOne <- function(vec) {
  p <- round(1024*vec/sum(vec),0)
  p[1] <- 1024-(sum(p)-p[1])
  p <- p/1024
  return(p)
}

并在原始代码中替换这些行：

  p <- priors/sum(priors) # turn arbitrary threesome into numbers that sum to 1
  p[1] <- 1-(sum(p)-p[1]) # ensure that numbers sum to 1

有了这个：

  p <- makeSumToOne(priors)

我正在阅读有关精度的文章，并注意到经常出现“2 的幂”，所以我认为在2^10=1024我的程序中使用“2 的幂”可能会有所帮助......到目前为止它很有帮助，但我怀疑它是否强大。除非没有人想出更好的答案，否则我不会将此视为答案。（或解释为什么此解决方案有效，以及具有round,signif或其他功能的解决方案失败。）

于 2013-11-23T02:02:50.910 回答