Проблема определения объема при использовании пакета LogisitcDx - PullRequest
Новая
       21

Проблема определения объема при использовании пакета LogisitcDx

0 голосов
/ 30 марта 2020

Я использую R, чтобы проверить, подходит ли модель, я также использую пакет paramtest для вычисления мощности и пакет LogisitcDx, чтобы оценить, подходит ли модель с помощью теста Штукеля. странная ошибка определения объема: 

power_model <- function(sim_num, N, m, p, groups=10, alpha=0.05){
    check_significance <- function(model, verbose=FALSE, groups=10, alpha=0.05){
        tests <- suppressWarnings(suppressMessages(LogisticDx::gof(model, g=groups, plotROC=FALSE)))
        tests <- tests[4]$gof
        ps <- tests$pVal
        p_HL <- ps[1]
        sig_HL <- (p_HL <= alpha)
        p_stukel <- ps[9]
        sig_stukel <- (p_stukel <= alpha)
        if(verbose){
            print(tests)
        }
        return(c(sig_HL=as.integer(sig_HL), sig_stukel=as.integer(sig_stukel)))
    }
    rate <- stats::rpois(1, lambda = N)
    probs <- rep(1 / m, m)
    xs <- stats::rmultinom(n=1, size=rate, prob=probs)
    xs <- xs[, 1]
    pis <- 1 - (1 - p) ** as.numeric(xs)
    ys <- Rlab::rbern(length(pis), pis)
    dd = data.frame(xs=xs, ys=ys)
    mo <- base::suppressWarnings(base::suppressMessages(stats::glm(ys ~ xs, data=dd, family=binomial('logit'), na.action = na.omit, maxit = 1000)))
    sigs <- base::suppressWarnings(base::suppressMessages(check_significance(mo)))
    sig_HL <- sigs['sig_HL']
    sig_stukel <- sigs['sig_stukel']
    return(c(exp_sig_HL=sig_HL, exp_sig_stukel=sig_stukel))

}
Ns <- c(300, 400, 500)
ms <- c(30, 50, 70, 100)
ps <- c(0.2, 0.4, 0.6, 0.8)
power_test <- grid_search(power_model, params=list(N=Ns, m=ms, p=ps), n.iter=100, output='data.frame',
                           # , parallel='snow', ncpus=12
                         )
power_result <- results(power_test) 
temp <- power_result %>%
    group_by(m.test, p.test) %>%
    summarise(
        power_HL=mean(exp_sig_HL, na.rm=TRUE),
        na_HL=sum(is.na(exp_sig_HL)),
        power_Stukel=mean(exp_sig_stukel, na.rm=TRUE),
        na_stukel=sum(is.na(exp_sig_stukel))
    )
temp

Обратите внимание, что я включил функцию 'check_significance' просто для ясности. Код отлично работает, когда я комментирую следующие строки: 

sigs <- base::suppressWarnings(base::suppressMessages(check_significance(mo)))
sig_HL <- sigs['sig_HL']
sig_stukel <- sigs['sig_stukel']

, но в противном случае вывод: error in is.data.frame(data): object 'dd' not found.

Я подозреваю, что это может быть связано с используемым пакетом (LogisticDx), я пытался исследовать это дальше, но я не могу понять, как использование пакета каким-то образом делает некоторые переменные несуществующими, и, чтобы добавить оскорбление ране, печать объекта 'dd' работает так, что это не так как-то удалено. Каковы ваши предложения?

...