purrr :: map () не работает хорошо с ordinal :: nominal_test

Question

Я пытался запустить номинальный тест для нескольких упорядоченных моделей логистической регрессии в столбце списка, но я получил пропущенные значения в выходных данных. Мне удалось запустить номинальный тест для отдельной модели, поместив подкадр данных. Я также пробовал циклы lapply () и for (), однако оба из них дали одинаковые результаты из purrr: map (). Я подозреваю, что это может быть связано с самой функцией nominal_test. Любые мысли по устранению проблемы приветствуются.

Ниже приведены мои воспроизводимые данные и коды.

Yanxian

library(tidyverse) 
library(ordinal) # for running ordered logistic regression
#> 
#> Attaching package: 'ordinal'
#> The following object is masked from 'package:dplyr':
#> 
#>     slice

# make a dataframe
df <- data.frame(stringsAsFactors=FALSE,
              Diet = c("REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM",
                       "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM",
                       "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "REF", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM",
                       "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM",
                       "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "REF", "REF", "REF", "REF", "REF", "REF", "REF",
                       "REF", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM",
                       "IM", "IM", "IM", "IM", "IM", "IM", "IM", "IM"),
           Segment = c("PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI",
                       "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI",
                       "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI", "PI",
                       "PI", "PI", "PI", "PI", "PI", "PI", "PI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI",
                       "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI", "MI"),
          Variable = c("hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv", "hpv",
                       "hpv", "hpv", "hpv", "hpv", "hpv", "lpc", "lpc", "lpc", "lpc",
                       "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc",
                       "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc",
                       "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc", "lpc",
                       "lpc", "lpc", "lpc", "lpc", "lpc", "lpc"),
              Rank = c("Severe", "Marked", "Moderate", "Moderate", "Moderate",
                       "Marked", "Marked", "Marked", "Moderate", "Severe",
                       "Moderate", "Marked", "Moderate", "Marked", "Marked", "Mild",
                       "Moderate", "Marked", "Moderate", "Moderate", "Moderate",
                       "Moderate", "Mild", "Moderate", "Mild", "Mild",
                       "Moderate", "Moderate", "Marked", "Moderate", "Moderate",
                       "Moderate", "Marked", "Moderate", "Moderate", "Marked", "Normal",
                       "Mild", "Mild", "Moderate", "Normal", "Normal", "Normal",
                       "Normal", "Mild", "Mild", "Normal", "Moderate", "Normal",
                       "Normal", "Normal", "Normal", "Mild", "Mild", "Normal",
                       "Normal", "Mild", "Normal", "Normal", "Mild", "Normal",
                       "Normal", "Normal", "Normal", "Normal", "Normal", "Normal",
                       "Normal", "Normal", "Normal", "Mild", "Normal", "Normal",
                       "Mild", "Normal", "Normal", "Normal", "Mild", "Moderate",
                       "Normal", "Normal", "Mild", "Normal", "Mild", "Normal",
                       "Normal", "Normal", "Normal", "Mild", "Normal", "Mild",
                       "Moderate", "Normal", "Normal", "Normal", "Normal", "Normal",
                       "Normal", "Normal", "Normal", "Normal", "Mild", "Normal",
                       "Normal", "Normal", "Mild")
      )

# Convert character to factor
df$Diet <- factor(df$Diet, levels = c("REF", "IM"))
df$Segment <- factor(df$Segment, levels = c("PI", "MI", "DI"))
df$Rank <- ordered(df$Rank, levels = c("Normal", "Mild", "Moderate", "Marked", "Severe")) # for running ordered logistic regression
df$Variable <- factor(df$Variable, levels = c("hpv", "snv", "smc", "lpc", "mfh"))

# run ordered logistic regression (olr), extract p values and check model assumption
olr <- df %>%
  group_by(Segment, Variable) %>% 
  nest() %>% 
  mutate(mod      = map(data, ~clm(Rank ~ Diet, data = .x)), # run olr for each dataframe
         Nominal  = map(mod, nominal_test), # somehow, nominal_test does NOT work
         p_nmn    = map_dbl(Nominal, ~.x["Diet", "Pr(>Chi)"]) # extact p value of nominal test
  )                                        

# nominal test can be done by subsetting dataframe directly
# example1
mod1 <- clm(Rank ~ Diet, data = df[1:36, ]) # Segment = "PI", Variable = "hpv"     
Nominal1 <- nominal_test(mod1)
p_nmn1 <- Nominal1["Diet", "Pr(>Chi)"]
p_nmn1
#> [1] 0.6356792

# example2
mod2 <- clm(Rank ~ Diet, data = df[73:107, ]) # Segment = "MI", Variable = "lpc"     
Nominal2 <- nominal_test(mod2)
p_nmn2 <- Nominal2["Diet", "Pr(>Chi)"]
p_nmn2
#> [1] 0.6742354

^{Создано в 2018-11-10 пакетом Представления (v0.2.1)}