मैंने प्रत्येक सेल के लिए परीक्षण आंकड़ों की गणना करने के लिए फॉर-लूप का उपयोग करके आर में एक ची-स्क्वायर टेस्ट लागू किया। हालाँकि, मैं सोच रहा था कि क्या इसे अनुकूलित किया जा सकता है। और क्या आर में ची = वर्ग मेरे कोड के रूप में काम कर रहा है?

eval_preds <- function(df, observations, predictions, groups) {
  # determine cells
  cells <- unique(df[,groups])
  my_sum = 0

  # compute test statistic for every cell
  for (i in 1:length(cells)) {

    # get cell means
    m_obs <- mean(df[df[,groups] == cells[i], observations])
    m_pred <- mean(df[df[,groups] == cells[i], predictions])

    # get cell variance
    var_obs <- var(df[df[,groups] == cells[i], observations])
    var_pred <- var(df[df[,groups] == cells[i], predictions])

    # get cell's number of observations
    N_obs <- length(df[df[,groups] == cells[i], observations])
    N_pred <- length(df[df[,groups] == cells[i], predictions])

    # sum up
    my_sum = my_sum + (m_obs-m_pred) / ((var_obs/N_obs) + (var_pred/N_pred))
  }
  return(my_sum)
}

और एक खिलौना उदाहरण:

# save this as df_head.csv
"RT","RTmodel_s","groups"
899,975.978308710825,"pl_sgDom"
1115,1095.61538562629,"sg_sgDom"
1077,1158.19266217845,"pl_sgDom"
850,996.410033287249,"sg_plDom"
854,894.862587823602,"pl_sgDom"
1720,1046.34200941684,"sg_sgDom"

# load dat into R
df <- read.csv('./df_head.csv')

my_chi <- eval_preds(df, 'RT', 'Pred', 'Group') 

संपादित करें

फ़ंक्शन eval_preds फ़ंक्शन को लूप के लिए कहा जाता है जिसमें एक मुक्त पैरामीटर t_parsing के आधार पर विभिन्न भविष्यवाणियों की गणना की जाती है।

p_grid = seq(0,1,0.1)

# tune p
for (i in 1:length(p_grid)) {
  
  # set t_parsing
  t_parsing = p_grid[i]
  
  # compute model-time RTs
  df$RTmodel <- ifelse(df$number == 'sing', 
                             RT_lookup(df$sg_t_act, df$epsilon), 
                             RT_decompose(df$sg_t_act, df$pl_t_act, df$epsilon))
  
  # scale into real time
  df$RTmodel_s <- scale_RTmodel(df$RT, df$RTmodel)
  
  # compare model output to measured RT
  # my_chi <- eval_preds(my_nouns, 'RT', 'RTmodel_s', 'groups')
  my_chi <- eval_preds1(df, RT, RTmodel_s, groups) #function written by DaveArmstrong
  print(paste(p_grid[i], ': ', my_chi))
}

RT और RTmodel_s संख्यात्मक चर हैं, groups एक वर्ण चर है।

-2
hyhno01 11 सितंबर 2020, 17:59

2 जवाब

सबसे बढ़िया उत्तर

ज़रूर, आप इसे dplyr और rlang के साथ कर सकते हैं:

eval_preds <- function(df, observations, predictions, groups) {
  # determine cells
  require(dplyr)
  require(rlang)
  groups <- enquo(groups)
  observations <- enquo(observations)
  predictions <- enquo(predictions)
  out <- DF %>% 
    group_by(!!groups) %>% 
    summarise(m_obs = mean(!!observations), 
              m_pred = mean(!!predictions), 
              var_obs = var(!!observations), 
              var_pred = var(!!predictions),
              N_obs = sum(!is.na(!!observations)),
              N_pred = sum(!is.na(!!predictions))) %>% 
    mutate(my_sum = (m_obs - m_pred)/((var_obs/N_obs) + (var_pred/N_pred))) 
sum(out$my_sum)
}
my_chi <- eval_preds(DF, RT, Pred, Group) 
my_chi
# [1] 1.6

या, थोड़ा और सुव्यवस्थित:

eval_preds <- function(df, observations, predictions, groups) {
  # determine cells
  require(dplyr)
  require(rlang)
  groups <- enquo(groups)
  observations <- enquo(observations)
  predictions <- enquo(predictions)
  out <- DF %>% 
    group_by(!!groups) %>% 
    summarise(diff= mean(!!observations) - mean(!!predictions), 
              sum_v = (var(!!observations)/n()) + (var(!!predictions)/n())) %>% 
    mutate(my_sum = diff/sum_v) 
  sum(out$my_sum)
}
my_chi <- eval_preds(DF, RT, Pred, Group) 
my_chi
# [1] 1.6

संपादित करें - अतिरिक्त मानदंड

इसलिए, मुझे लगता है कि कौन सा प्रश्न बेहतर है यह डेटा सेट के आकार पर निर्भर करता है। मैं मिश्रण में एक और फ़ंक्शन भी फेंकना चाहता हूं - वह जो आधार आर से by() का उपयोग करता है:

eval_preds <- function(df, observations, predictions, groups) {
  # determine cells
  cells <- unique(df[,groups])
  my_sum = 0
  
  # compute test statistic for every cell
  for (i in 1:length(cells)) {
    
    # get cell means
    m_obs <- mean(df[df[,groups] == cells[i], observations])
    m_pred <- mean(df[df[,groups] == cells[i], predictions])
    
    # get cell variance
    var_obs <- var(df[df[,groups] == cells[i], observations])
    var_pred <- var(df[df[,groups] == cells[i], predictions])
    
    # get cell's number of observations
    N_obs <- length(df[df[,groups] == cells[i], observations])
    N_pred <- length(df[df[,groups] == cells[i], predictions])
    
    # sum up
    my_sum = my_sum + (m_obs-m_pred) / ((var_obs/N_obs) + (var_pred/N_pred))
  }
  return(my_sum)
}

eval_preds1 <- function(df, observations, predictions, groups) {
  # determine cells
  require(dplyr)
  require(rlang)
  groups <- enquo(groups)
  observations <- enquo(observations)
  predictions <- enquo(predictions)
  out <- df %>% 
    group_by(!!groups) %>% 
    summarise(m_obs = mean(!!observations), 
              m_pred = mean(!!predictions), 
              var_obs = var(!!observations), 
              var_pred = var(!!predictions),
              N_obs = sum(!is.na(!!observations)),
              N_pred = sum(!is.na(!!predictions))) %>% 
    ungroup %>%
    mutate(my_sum = (m_obs - m_pred)/((var_obs/N_obs) + (var_pred/N_pred))) 
  sum(out$my_sum)
}

eval_preds2 <- function(df, observations, predictions, groups) {
  # determine cells
  require(dplyr)
  require(rlang)
  groups <- enquo(groups)
  observations <- enquo(observations)
  predictions <- enquo(predictions)
  out <- df %>% 
    group_by(!!groups) %>% 
    summarise(diff= mean(!!observations) - mean(!!predictions), 
              sum_v = (var(!!observations)/n()) + (var(!!predictions)/n())) %>% 
    ungroup %>%
    mutate(my_sum = diff/sum_v) 
  sum(out$my_sum)
}

eval_preds3<- function(df, observations, predictions, groups) {
  # determine cells
  
  m <- by(df[,c(observations, predictions)], list(df[[groups]]), function(x)diff(-colMeans(x)))
  v <- by(df[,c(observations, predictions)], list(df[[groups]]), function(x)sum(apply(x, 2, var)/nrow(x)))
  sum(m/v)
}

तो, eval_preds() मूल है, eval_preds1() dplyr कोड का पहला सेट है और eval_preds2(), eval_preds3() by() के साथ आधार R है। मूल डेटा सेट पर, यहाँ microbenchmark() से आउटपुट है।

microbenchmark(eval_preds(DF, 'RT', 'Pred', 'Group'), 
               eval_preds1(DF, RT, Pred, Group), 
               eval_preds2(DF, RT, Pred, Group), 
               eval_preds4(DF, 'RT', 'Pred', 'Group'), times=100)
# Unit: microseconds
#                                  expr      min        lq      mean   median        uq       max neval  cld
#  eval_preds(DF, "RT", "Pred", "Group")  236.513  279.4920  324.9760  295.190  321.0125   774.213   100 a   
#       eval_preds1(DF, RT, Pred, Group) 5236.850 5747.5095 6503.0251 6089.343 6937.8670 12950.677   100    d
#       eval_preds2(DF, RT, Pred, Group) 4871.812 5372.2365 6070.7297 5697.686 6548.8935 14577.786   100   c 
# eval_preds4(DF, "RT", "Pred", "Group")  651.013  739.9405  839.1706  773.610  923.9870  1582.218   100  b  

इस मामले में, मूल कोड सबसे तेज़ है। हालांकि, अगर हमने एक बड़ा डेटासेट बनाया है - यहां 1000 समूहों वाला एक है और प्रति समूह 10 अवलोकन हैं।

DF2 <- data.frame(
  Group = rep(1:1000, each=10), 
  RT = rpois(10000, 3), 
  Pred = rpois(10000, 4)
)

यहाँ से microbenchmark() का आउटपुट बिल्कुल अलग कहानी कहता है।

microbenchmark(eval_preds(DF2, 'RT', 'Pred', 'Group'), 
               eval_preds1(DF2, RT, Pred, Group), 
               eval_preds2(DF2, RT, Pred, Group), 
               eval_preds3(DF2, 'RT', 'Pred', 'Group'), times=25)


# Unit: milliseconds
#                                   expr       min        lq     mean    median        uq      max neval cld
#  eval_preds(DF2, "RT", "Pred", "Group") 245.67280 267.96445 353.6489 324.29416 419.4193 565.4494    25   c
#       eval_preds1(DF2, RT, Pred, Group)  74.56522  88.15003 102.6583  92.24103 106.6766 211.2368    25 a  
#       eval_preds2(DF2, RT, Pred, Group)  79.00919  89.03754 125.8202  94.71703 114.5176 606.8830    25 a  
# eval_preds3(DF2, "RT", "Pred", "Group") 193.94042 240.35447 272.0004 254.85557 316.5098 420.5394    25  b 

दोनों dplyr() आधारित फंक्शन अपने बेस R प्रतिस्पर्धियों की तुलना में बहुत तेज हैं। तो, "इष्टतम" पद्धति का प्रश्न हल की जाने वाली समस्या के आकार पर निर्भर करता है।

2
DaveArmstrong 13 सितंबर 2020, 17:46

यहां एक आधार दृष्टिकोण है जो आपके लूप के समान दिखना चाहिए। ध्यान दें, लूप अक्षम है क्योंकि प्रत्येक लूप फ़िल्टरिंग के लिए प्रत्येक समूह को देख रहा है। ग्रुपिंग ऑपरेशंस, जैसे by() मदद करेंगे क्योंकि हमें ग्रुपिंग प्राप्त करने के लिए ग्रुपिंग वेरिएबल को केवल एक बार स्कैन करना होगा।

DF <- read.table(header=T, text="Group    Pred    RT
cond1       2        3
cond1       4        2
cond2       2        2
cond2       1        2")

stats = by(data = DF[c("Pred", "RT")], 
           INDICES = DF["Group"],
           simplify = FALSE,
           FUN = function(DF_grp) {
             RT = DF_grp[["RT"]]
             Pred = DF_grp[["Pred"]]
             
             m_obs = mean(RT)
             m_pred = mean(Pred)
             
             var_obs = var(RT)
             var_pred = var(Pred)
             
             N_obs = N_pred = length(Pred) ##simplification
             
             return((m_obs - m_pred) / ((var_obs / N_obs) + (var_pred / N_pred)))
           }
)
           
stats
#> Group: cond1
#> [1] -0.4
#> ------------------------------------------------------------ 
#> Group: cond2
#> [1] 2

sum(unlist(stats, use.names = FALSE))
#> [1] 1.6

अंत में, यदि आप जानते हैं कि आपके समूहीकरण का आदेश दिया गया है, तो आप अतिरिक्त प्रदर्शन के लिए Rcpp देख सकते हैं।

0
Cole 13 सितंबर 2020, 16:19