带有函数错误的 R 和 nls -“必须具有相同的行数”答案

【问题标题】：R and nls with functions error - 'must have same number of rows'带有函数错误的 R 和 nls -“必须具有相同的行数”
【发布时间】：2014-05-31 22:16:43
【问题描述】：

我想为从危险率（死亡率）确定生存的函数找到最佳参数。

eh <- function(kappa,lambda,time,delay){
  FoM <- -log(1-(1-exp(-kappa*(time+delay)))*(exp(-lambda*time)))  
  return(FoM)  
}

survival <- function(k, l, t, d){
  theSurvs <- array(1, c(1, length(t)))
  S = array(1)
  Qs <- array(0)
  for(j in 1:length(t)){
    S[1] = 1
    for(i in 1:(t[j]+1)){
      theEHs <- eh(k, l, i-1, d)
      theQs <- hazards2Qs(theEHs)
      S[i+1] <- S[i]-S[i]*theQs
    }
    theSurvs[j] <- S[i]
  }
  return(theSurvs)
}

hazards2Qs <- function(hazards){
  Qs<- 1-exp(-hazards)
  return(Qs)
}

timeX = 0  # time cycles of arbitrary length
kappaX = 0.001
lambdaX = 0.05
delayX = 50

theYears <- floor(runif(20)*100)
theYears

EHs <- eh(kappaX, lambdaX, theYears, delayX)


theS <- survival(kappaX, lambdaX, theYears, delayX)
theS

theData <- data.frame(theYears, t(theS))
theData
survival(kappaX, lambdaX, theYears, delayX)

mod <- nls(~ survival(kappa, lambda, theYears, delay), start=list(lambda=0.0015, kappa=0.0013, delay=48), data=theData, trace=1)

当我运行它时，我得到了这个错误：

3.647752 :   0.0015  0.0013 48.0000
Error in qr.qty(QR, resid) : 
  'qr' and 'y' must have the same number of rows

回溯（）：

4: stop("'qr' and 'y' must have the same number of rows")
3: qr.qty(QR, resid)
2: (function () 
   {
       if (npar == 0) 
           return(0)
       rr <- qr.qty(QR, resid)
       sqrt(sum(rr[1L:npar]^2)/sum(rr[-(1L:npar)]^2))
   })()

> dim(theS)
[1]  1 20
> dim(theYears)
NULL
> dim(theData)
[1] 20  2
> typeof(theData)
[1] "list"
> typeof(theS)
[1] "double"
> typeof(theYears)
[1] "double"

我已经苦苦挣扎了一天，并没有深究。有什么想法吗？

【问题讨论】：

没有样本数据，这个错误不是reproducible。如果您有一个返回相同错误的完整运行示例，那么获得帮助会容易得多。有关包含示例数据集的建议，请参阅链接。
抱歉，我未能正确复制其中一个函数。现在放回去，它应该会运行以重现提到的错误。
样本数据由代码生成为进入数据框“theData”的虚拟数据。
R 版本 3.1.0，在 Ubuntu 12.04 中运行并使用 RStudio。没有添加任何包。

标签： r nls

【解决方案1】：

在花了很长时间试图弄清楚自己发生了什么之后，我认为问题在于您如何塑造您的 survival 函数结果。我不确定你为什么要返回一个数组，但你应该在这里返回一个向量。所以只需将第一行更改为

#OLD: theSurvs <- array(1, c(1, length(t)))
theSurvs <- rep.int(1, length(t)

这意味着你也必须改变

#OLD: theData <- data.frame(theYears, t(theS))
theData <- data.frame(theYears, theS)

通过返回这样的形状对象，它会干扰使用qr() 的梯度计算。尽量远离一维数组，只使用简单的向量。

现在，话虽如此，解决这个问题似乎会导致另一个问题。似乎当它试图取数值导数时，它遇到了NaN/Inf 值。您可以在return(theSurvs) 行上方添加此代码，以查看发生这种情况时调用的参数

if(any(!is.finite(theSurvs))) {
    dput(c(k,l,d))
    dput(t)
    print(theSurvs)     
}

【讨论】：

太棒了。谢谢。我知道这将是一件简单而愚蠢的事情。