Error in eval(object$data)[, vars] : object of type 'closure' is not subsettable

Trying to generate confidence intervals for a logistic function using a bootstrap procedure, and have now spent several hours trying to figure out what is causing this error. Sometimes it works, but then it will fail, even on the same data. I am totally mystified. Reprex below:


library(tidyverse)
library(lubridate)
library(car)

indata <- tribble(
~Date,       ~Cases,  ~Days, 
"2020-03-10",    21,     0, 
"2020-03-11",    23,     1, 
"2020-03-12",    38,     2, 
"2020-03-13",    51,     3, 
"2020-03-14",    56,     4, 
"2020-03-15",    57,     5, 
"2020-03-16",    64,     6, 
"2020-03-17",    83,     7, 
"2020-03-18",   143,     8, 
"2020-03-19",   194,     9, 
"2020-03-20",   304,    10, 
"2020-03-21",   334,    11, 
"2020-03-22",   352,    12, 
"2020-03-23",   410,    13, 
"2020-03-24",   974,    14, 
"2020-03-25",  1396,    15, 
"2020-03-26",  1731,    16, 
"2020-03-27",  2052,    17, 
"2020-03-28",  2371,    18, 
"2020-03-29",  2877,    19, 
"2020-03-30",  3266,    20, 
"2020-03-31",  3997,    21, 
"2020-04-01",  4669,    22, 
"2020-04-02",  5330,    23, 
"2020-04-03",  6110,    24, 
"2020-04-04",  6812,    25, 
"2020-04-05",  7276,    26 
)

indata$Date <- as_date(indata$Date)


#---------------------------------------------------    
#------------------- Fit Logistic function ---------
#---------------------------------------------------    
logistic_model <- function(df, 
                     indep="Cases", # independent variable
                     r=0.24,
                     projection=10){
    
    print(df)
    Asym <- max(df$Cases)*5
    xmid <- max(df$Days)*2
    scal <- 1/r
    
    ## using a selfStart model
    my_model   <- nls(Cases ~ SSlogis(Days, Asym, xmid, scal), 
                   data=df)
    coeffs <- coef(my_model)
    
    dayseq <- df$Days
    dayseq <- c(dayseq,(dayseq[length(dayseq)]+1):
                    (dayseq[length(dayseq)]-1+projection))
    dateseq <- df$Date
    dateseq <- as_date(c(dateseq,(dateseq[length(dateseq)]+1): 
                             (dateseq[length(dateseq)]-1+projection)))
    foo <- tibble(Days=dayseq)
    predict2 <- function(x) {stats::predict(x, foo)}
    print(summary(foo))
    f.boot <- car::Boot(my_model,f=predict2)
    
    Cases <- predict(my_model, data.frame(Days=dayseq))
    confidence <- confint(f.boot)
    preds <- tibble(dayseq, dateseq,
                    Cases,
                    confidence["2.5 %"],
                    confidence["97.5 %"])
    names(preds) <- c("Days", "Date","Cases","SD_upper","SD_lower")
    
    #  return a tibble
    tribble(~Line, ~r, ~K, ~xmid,
            preds, 1/coeffs[["scal"]], coeffs[["Asym"]], coeffs[["xmid"]])
    
} 

        foo <- logistic_model(indata, 
                        indep="Cases", # independent variable
                        r=0.24,
                        projection=10)

This is exactly what happens whenever an object is a function, f is asked to return a subset

mean[1]
#> Error in mean[1]: object of type 'closure' is not subsettable

Created on 2020-04-06 by the reprex package (v0.3.0)

Once the head-fake term closure is translated to function and taking operator precedence into account

f(obj) returns a value which, itself may or may not be subsetable

Here it is

a <- mean(seq(1:10))
a[1]
#> [1] 5.5

Created on 2020-04-06 by the reprex package (v0.3.0)

Here, the subset operator [\ ] is bound to a, after mean has done its work. Until that, it's bound to mean(a) and there's nothing that can be subset.

It appears to be an issue with argument handling in the logistic_model function. Unpacking the pieces seems to run as desired.

# as in OP
library(tidyverse)
library(lubridate)
#> 
#> Attaching package: 'lubridate'
#> The following objects are masked from 'package:dplyr':
#> 
#>     intersect, setdiff, union
#> The following objects are masked from 'package:base':
#> 
#>     date, intersect, setdiff, union
library(car)
#> Loading required package: carData
#> 
#> Attaching package: 'car'
#> The following object is masked from 'package:dplyr':
#> 
#>     recode
#> The following object is masked from 'package:purrr':
#> 
#>     some

indata <- tribble(
  ~Date,       ~Cases,  ~Days, 
  "2020-03-10",    21,     0, 
  "2020-03-11",    23,     1, 
  "2020-03-12",    38,     2, 
  "2020-03-13",    51,     3, 
  "2020-03-14",    56,     4, 
  "2020-03-15",    57,     5, 
  "2020-03-16",    64,     6, 
  "2020-03-17",    83,     7, 
  "2020-03-18",   143,     8, 
  "2020-03-19",   194,     9, 
  "2020-03-20",   304,    10, 
  "2020-03-21",   334,    11, 
  "2020-03-22",   352,    12, 
  "2020-03-23",   410,    13, 
  "2020-03-24",   974,    14, 
  "2020-03-25",  1396,    15, 
  "2020-03-26",  1731,    16, 
  "2020-03-27",  2052,    17, 
  "2020-03-28",  2371,    18, 
  "2020-03-29",  2877,    19, 
  "2020-03-30",  3266,    20, 
  "2020-03-31",  3997,    21, 
  "2020-04-01",  4669,    22, 
  "2020-04-02",  5330,    23, 
  "2020-04-03",  6110,    24, 
  "2020-04-04",  6812,    25, 
  "2020-04-05",  7276,    26 
)

indata$Date <- as_date(indata$Date)

# Unpack function to inspect steps
# logistic_model <- function(df., 
#                           indep="Cases", # independent variable
#                           r=0.24,
#                           projection=10){

df. <- indata           # avoid potential namespace confict with df.()
df.
#> # A tibble: 27 x 3
#>    Date       Cases  Days
#>    <date>     <dbl> <dbl>
#>  1 2020-03-10    21     0
#>  2 2020-03-11    23     1
#>  3 2020-03-12    38     2
#>  4 2020-03-13    51     3
#>  5 2020-03-14    56     4
#>  6 2020-03-15    57     5
#>  7 2020-03-16    64     6
#>  8 2020-03-17    83     7
#>  9 2020-03-18   143     8
#> 10 2020-03-19   194     9
#> # … with 17 more rows
r <- 0.24
projection <- 10
indep <- df.$Cases
Asym <- max(indep)
Asym
#> [1] 7276
Asym <- max(indep)*5
Asym
#> [1] 36380
xmid <- max(indep)*2
xmid
#> [1] 14552
scal <- 1/r
scal
#> [1] 4.166667

## using a selfStart model
my_model   <- nls(Cases ~ SSlogis(Days, Asym, xmid, scal), data=df.)
my_model
#> Nonlinear regression model
#>   model: Cases ~ SSlogis(Days, Asym, xmid, scal)
#>    data: df.
#>      Asym      xmid      scal 
#> 10130.468    22.497     3.706 
#>  residual sum-of-squares: 357740
#> 
#> Number of iterations to convergence: 0 
#> Achieved convergence tolerance: 1.495e-06

coeffs <- coef(my_model)
coeffs
#>         Asym         xmid         scal 
#> 10130.467944    22.496894     3.705809

dayseq <- df.$Days
dayseq
#>  [1]  0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
#> [26] 25 26
dayseq <- c(dayseq,(dayseq[length(dayseq)]+1):
              (dayseq[length(dayseq)]-1+projection))
dayseq
#>  [1]  0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
#> [26] 25 26 27 28 29 30 31 32 33 34 35
dateseq <- df.$Date
dateseq
#>  [1] "2020-03-10" "2020-03-11" "2020-03-12" "2020-03-13" "2020-03-14"
#>  [6] "2020-03-15" "2020-03-16" "2020-03-17" "2020-03-18" "2020-03-19"
#> [11] "2020-03-20" "2020-03-21" "2020-03-22" "2020-03-23" "2020-03-24"
#> [16] "2020-03-25" "2020-03-26" "2020-03-27" "2020-03-28" "2020-03-29"
#> [21] "2020-03-30" "2020-03-31" "2020-04-01" "2020-04-02" "2020-04-03"
#> [26] "2020-04-04" "2020-04-05"
dateseq <- as_date(c(dateseq,(dateseq[length(dateseq)]+1): 
                       (dateseq[length(dateseq)]-1+projection)))
dateseq
#>  [1] "2020-03-10" "2020-03-11" "2020-03-12" "2020-03-13" "2020-03-14"
#>  [6] "2020-03-15" "2020-03-16" "2020-03-17" "2020-03-18" "2020-03-19"
#> [11] "2020-03-20" "2020-03-21" "2020-03-22" "2020-03-23" "2020-03-24"
#> [16] "2020-03-25" "2020-03-26" "2020-03-27" "2020-03-28" "2020-03-29"
#> [21] "2020-03-30" "2020-03-31" "2020-04-01" "2020-04-02" "2020-04-03"
#> [26] "2020-04-04" "2020-04-05" "2020-04-06" "2020-04-07" "2020-04-08"
#> [31] "2020-04-09" "2020-04-10" "2020-04-11" "2020-04-12" "2020-04-13"
#> [36] "2020-04-14"
foo <- tibble(Days=dayseq)
foo
#> # A tibble: 36 x 1
#>     Days
#>    <dbl>
#>  1     0
#>  2     1
#>  3     2
#>  4     3
#>  5     4
#>  6     5
#>  7     6
#>  8     7
#>  9     8
#> 10     9
#> # … with 26 more rows
predict2 <- function(x) {stats::predict(x, foo)}
print(summary(foo))
#>       Days      
#>  Min.   : 0.00  
#>  1st Qu.: 8.75  
#>  Median :17.50  
#>  Mean   :17.50  
#>  3rd Qu.:26.25  
#>  Max.   :35.00
f.boot <- car::Boot(my_model,f=predict2)
#> Loading required namespace: boot
f.boot
#> 
#> ORDINARY NONPARAMETRIC BOOTSTRAP
#> 
#> 
#> Call:
#> boot::boot(data = data.frame(update(object, model = TRUE)$model), 
#>     statistic = boot.f, R = R, .fn = f, parallel = p_type, ncpus = ncores, 
#>     cl = cl2)
#> 
#> 
#> Bootstrap Statistics :
#>        original       bias    std. error
#> t1*    23.34271   -0.4866764    6.170080
#> t2*    30.55162   -0.6494386    7.634477
#> t3*    39.97806   -0.8572582    9.413758
#> t4*    52.29787   -1.1194215   11.561217
#> t5*    68.38849   -1.4456372   14.132285
#> t6*    89.38583   -1.8450955   17.180618
#> t7*   116.75519   -2.3248455   20.751736
#> t8*   152.37769   -2.8872076   24.873422
#> t9*   198.65316   -3.5259009   29.542128
#> t10*  258.61775   -4.2206262   34.705157
#> t11*  336.07050   -4.9301107   40.239620
#> t12*  435.69562   -5.5842901   45.931675
#> t13*  563.15564   -6.0776290   51.463497
#> t14*  725.11438   -6.2677814   56.419930
#> t15*  929.12766   -5.9866532   60.328484
#> t16* 1183.31907   -5.0732167   62.738629
#> t17* 1495.74780   -3.4361798   63.324491
#> t18* 1873.39421   -1.1456286   61.971226
#> t19* 2320.76197    1.4664003   58.829197
#> t20* 2838.23972    3.7410084   54.471486
#> t21* 3420.56561    4.7221588   50.586974
#> t22* 4055.91561    3.3230391   51.475973
#> t23* 4726.15538   -1.3874649   63.485987
#> t24* 5408.53856   -9.8738393   89.480732
#> t25* 6078.62835  -21.9394712  127.932482
#> t26* 6713.69899  -36.7506952  176.224488
#> t27* 7295.65160  -53.0699325  231.491267
#> t28* 7812.70357  -69.5831710  290.621872
#> t29* 8259.63192  -85.1818215  350.520681
#> t30* 8636.85592  -99.1149886  408.479233
#> t31* 8948.89958 -111.0078632  462.444355
#> t32* 9202.75391 -120.7925148  511.108201
#> t33* 9406.48122 -128.6062318  553.848349
#> t34* 9568.20315 -134.6962370  590.584056
#> t35* 9695.47071 -139.3483623  621.609453
#> t36* 9794.94166 -142.8420482  647.442477
Cases <- predict(my_model, data.frame(Days=dayseq))
Cases
#>  [1]   23.34271   30.55162   39.97806   52.29787   68.38849   89.38583
#>  [7]  116.75519  152.37769  198.65316  258.61775  336.07050  435.69562
#> [13]  563.15564  725.11438  929.12766 1183.31907 1495.74780 1873.39421
#> [19] 2320.76197 2838.23972 3420.56561 4055.91561 4726.15538 5408.53856
#> [25] 6078.62835 6713.69899 7295.65160 7812.70357 8259.63192 8636.85592
#> [31] 8948.89958 9202.75391 9406.48122 9568.20315 9695.47071 9794.94166
#> attr(,"gradient")
#>              Asym        xmid       scal
#>  [1,] 0.002304208   -6.284438   38.15101
#>  [2,] 0.003015816   -8.219390   47.67957
#>  [3,] 0.003946319  -10.745372   59.43284
#>  [4,] 0.005162434  -14.039551   73.86448
#>  [5,] 0.006750773  -18.329821   91.49009
#>  [6,] 0.008823465  -23.907639  112.87939
#>  [7,] 0.011525152  -31.142883  138.63664
#>  [8,] 0.015041525  -40.500121  169.36279
#>  [9,] 0.019609475  -52.554703  205.59075
#> [10,] 0.025528707  -68.005554  247.68245
#> [11,] 0.033174233  -87.679007  295.67508
#> [12,] 0.043008440 -112.514452  349.06464
#> [13,] 0.055590289 -143.517840  406.52168
#> [14,] 0.071577580 -181.664114  465.55152
#> [15,] 0.091716164 -227.726715  522.14509
#> [16,] 0.116807938 -282.016183  570.52203
#> [17,] 0.147648441 -344.028278  603.13834
#> [18,] 0.184926720 -412.043283  611.19138
#> [19,] 0.229087342 -482.783908  585.84463
#> [20,] 0.280168669 -551.311217  520.23111
#> [21,] 0.337651294 -611.366485  411.92558
#> [22,] 0.400368041 -656.282299  265.09332
#> [23,] 0.466528833 -680.355590   91.22562
#> [24,] 0.533888325 -680.278775  -92.35562
#> [25,] 0.600034311 -656.062687 -266.10432
#> [26,] 0.662723482 -611.033443 -412.72542
#> [27,] 0.720169260 -550.904750 -520.77100
#> [28,] 0.771208558 -482.345390 -586.12103
#> [29,] 0.815325803 -411.608108 -611.23583
#> [30,] 0.852562386 -343.622014 -603.00209
#> [31,] 0.883364878 -281.654045 -570.26153
#> [32,] 0.908423378 -227.415172 -521.81197
#> [33,] 0.928533733 -181.403348 -465.18733
#> [34,] 0.944497648 -143.304155 -406.15662
#> [35,] 0.957060499 -112.342194 -348.71853
#> [36,] 0.966879488  -87.541887 -295.35942
confidence <- confint(f.boot)
confidence
#> Bootstrap bca confidence intervals
#> 
#>          2.5 %      97.5 %
#> V1    12.29421    36.78035
#> V2    16.62235    46.77471
#> V3    22.50012    59.62103
#> V4    30.35387    76.08839
#> V5    41.50583    97.16222
#> V6    56.45748   123.92827
#> V7    76.58277   157.73101
#> V8   103.19305   200.37641
#> V9   138.81629   254.28667
#> V10  186.76832   322.77210
#> V11  252.09683   409.49938
#> V12  339.23453   518.00073
#> V13  453.94899   655.36019
#> V14  601.14273   825.64160
#> V15  793.03875  1036.14563
#> V16 1038.62081  1288.26522
#> V17 1346.21539  1606.07691
#> V18 1725.10013  1976.35854
#> V19 2180.92446  2421.26546
#> V20 2716.17593  2935.47340
#> V21 3308.69274  3515.97235
#> V22 3974.49671  4167.46224
#> V23 4650.81530  4848.08154
#> V24 5297.29771  5530.97005
#> V25 5867.97492  6179.78124
#> V26 6421.29890  6796.79802
#> V27 7071.33555  7490.84967
#> V28 7477.23413  8129.73858
#> V29 7754.22767  8669.90743
#> V30 7955.41737  9126.02168
#> V31 8151.33134  9552.85961
#> V32 8284.31206  9913.10172
#> V33 8360.49797 10191.01474
#> V34 8422.28792 10417.73427
#> V35 8501.58178 10634.02729
#> V36 8551.54545 10792.10457
preds <- tibble(dayseq, dateseq,Cases, confidence["2.5 %"], confidence["97.5 %"])
preds
#> # A tibble: 36 x 5
#>    dayseq dateseq    Cases `confidence["2.5 %"]` `confidence["97.5 %"]`
#>     <dbl> <date>     <dbl>                 <dbl>                  <dbl>
#>  1      0 2020-03-10  23.3                    NA                     NA
#>  2      1 2020-03-11  30.6                    NA                     NA
#>  3      2 2020-03-12  40.0                    NA                     NA
#>  4      3 2020-03-13  52.3                    NA                     NA
#>  5      4 2020-03-14  68.4                    NA                     NA
#>  6      5 2020-03-15  89.4                    NA                     NA
#>  7      6 2020-03-16 117.                     NA                     NA
#>  8      7 2020-03-17 152.                     NA                     NA
#>  9      8 2020-03-18 199.                     NA                     NA
#> 10      9 2020-03-19 259.                     NA                     NA
#> # … with 26 more rows
names(preds) <- c("Days", "Date","Cases","SD_upper","SD_lower")
preds
#> # A tibble: 36 x 5
#>     Days Date       Cases SD_upper SD_lower
#>    <dbl> <date>     <dbl>    <dbl>    <dbl>
#>  1     0 2020-03-10  23.3       NA       NA
#>  2     1 2020-03-11  30.6       NA       NA
#>  3     2 2020-03-12  40.0       NA       NA
#>  4     3 2020-03-13  52.3       NA       NA
#>  5     4 2020-03-14  68.4       NA       NA
#>  6     5 2020-03-15  89.4       NA       NA
#>  7     6 2020-03-16 117.        NA       NA
#>  8     7 2020-03-17 152.        NA       NA
#>  9     8 2020-03-18 199.        NA       NA
#> 10     9 2020-03-19 259.        NA       NA
#> # … with 26 more rows

#  return a tibble
tribble(~Line, ~r, ~K, ~xmid,
        preds, 1/coeffs[["scal"]], coeffs[["Asym"]], coeffs[["xmid"]])
#> # A tibble: 1 x 4
#>   Line                  r      K  xmid
#>   <list>            <dbl>  <dbl> <dbl>
#> 1 <tibble [36 × 5]> 0.270 10130.  22.5

#} 

Created on 2020-04-06 by the reprex package (v0.3.0)

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.

That helps, at least I understand what it is objecting to. However, the error is coming from

f.boot <- car::Boot(my_model,f=predict2)

with a traceback of

Error in eval(object$data)[, vars] : object of type 'closure' is not subsettable
9. na.omit(eval(object$data)[, vars])
8. nls(formula = Cases ~ SSlogis(Days, Asym, xmid, scal), data = na.omit(eval(object$data)[, vars]), algorithm = "default", control = structure(list(maxiter = 50, tol = 1e-05, minFactor = 0.0009765625, printEval = FALSE, warnOnly = FALSE), .Names = c("maxiter", "tol", "minFactor", ...
7. eval(call, parent.frame())
6. eval(call, parent.frame())
5. update.default(object, data = na.omit(eval(object$data)[, vars]), start = coef(object))
4. update(object, data = na.omit(eval(object$data)[, vars]), start = coef(object))
3. Boot.nls(my_model, f = predict2)
2. car::Boot(my_model, f = predict2)

  1. logistic_model(indata, indep = "Cases", r = 0.24, projection = 10)

and I just don't see what it doesn't like about what I have fed it.

I have to go all Bill Clinton here: I feel your pain. I've got a similar hang-up about how to reliably pass data frames to functions and get them unpacked. I think I'm finally going to have to get initiated into the mysteries of quasiquotation and the other magic of tidyeval

You are correct there. And now I think I'm going crazy. I agree that works, but what could possibly be wrong with the function call? I've tried fiddling it in multiple ways to no avail.

1 Like

So it is something funky about passing a data frame as a function argument. Totally mysterious to me, so sadly I'll just have to make it a global for now, since that works. I was trying to be good, for my shiny app, and isolate data to each function, but c'est la vie. Truly odd how sometimes it works just fine for other functions I have in the same app. I'll mark this as the "solution", but not very satisfying. Thanks for your help! Sometimes R's eccentricities drive me to distraction.

1 Like