Построить эквивалент матрицы корреляции для факторов (категориальные данные)?А смешанные типы?

Question

На самом деле есть 2 вопроса, один из которых более сложный, чем другой.

В1: Я ищу метод, который похож на corrplot(), но может иметь дело с факторами.

Первоначально я пытался использовать chisq.test(), затем вычислил p-значение и V Крамера Крамера как корреляцию, но там слишком много столбцов, чтобы вычислить.Так может кто-нибудь сказать мне, если есть быстрый способ создать «corrplot», что каждая ячейка содержит значение Cramer's V , в то время как цвет отображается как p-значение .Или любой другой вид подобного сюжета.

Относительно V Cramer V , скажем, tbl - это двумерный факторный фрейм данных.

chi2 <- chisq.test(tbl, correct=F)
Cramer_V <- sqrt(chi2$/nrow(tbl)) 

Я подготовил тестовый фрейм данных с коэффициентами:

df <- data.frame(
group = c('A', 'A', 'A', 'A', 'A', 'B', 'C'),
student = c('01', '01', '01', '02', '02', '01', '02'),
exam_pass = c('Y', 'N', 'Y', 'N', 'Y', 'Y', 'N'),
subject = c('Math', 'Science', 'Japanese', 'Math', 'Science', 'Japanese', 'Math')
) 

Q2: Тогда я бы хотел вычислить матрицу корреляции / ассоциации на кадре данных смешанного типа, например:

df <- data.frame(
group = c('A', 'A', 'A', 'A', 'A', 'B', 'C'),
student = c('01', '01', '01', '02', '02', '01', '02'),
exam_pass = c('Y', 'N', 'Y', 'N', 'Y', 'Y', 'N'),
subject = c('Math', 'Science', 'Japanese', 'Math', 'Science', 'Japanese', 'Math')
) 
df$group <- factor(df$group, levels = c('A', 'B', 'C'), ordered = T)
df$student <- as.integer(df$student)

Answer 1

Что касается Q1, вы можете использовать ? Пар.таблицу из пакета vcd , если вы сначала преобразуете свой фрейм данных с ? Структурируемым (из того же самогопакет).Это даст вам матрицу графиков из мозаичных графиков .Это не совсем то, что делает corrplot(), но я подозреваю, что это будет более полезная визуализация.

df <- data.frame(
  ... 
) 
library(vcd)
st <- structable(~group+student+exam_pass+subject, df)
st
#                 student       01                    02             
#                 subject Japanese Math Science Japanese Math Science
# group exam_pass                                                    
# A     N                        0    0       1        0    1       0
#       Y                        1    1       0        0    0       1
# B     N                        0    0       0        0    0       0
#       Y                        1    0       0        0    0       0
# C     N                        0    0       0        0    1       0
#       Y                        0    0       0        0    0       0
pairs(st)

Существует множество других графиков, подходящих для категориально-категориальных данных, таких как графики сит, графики ассоциации и давлениеграфики (см. мой вопрос по Перекрестная проверка здесь: Альтернатива решетчатым / мозаичным графикам для таблиц сопряженности ).Вы можете написать свою собственную парную функцию, чтобы поместить все, что вы хотите, в верхнюю или нижнюю панель треугольника (см. Мой вопрос здесь: Матрица пар с qq-plots ), если вы не предпочитаете мозаичные графики.Просто помните, что хотя матрицы графиков очень полезны, они всегда отображают только маргинальные проекции (чтобы понять это более подробно, посмотрите мои ответы по резюме здесь: Есть ли разница между «контролем» и «игнорированием» других переменных в несколькихрегрессия? , а здесь: Альтернативы трехмерному точечному графику ).

Что касается Q2, вам нужно написать пользовательскую функцию.

Answer 2

Решение от @AntoniosK может быть улучшено в соответствии с предложением @JD, чтобы учесть смешанные кадры данных, включающие как номинальные , так и числовые атрибуты.Сила ассоциации рассчитывается для номинального и номинального с исправленным смещением V Крамера, числового и числового с корреляцией Спирмена (по умолчанию) или Пирсона, и номинального против числового с ANOVA.

require(tidyverse)
require(rcompanion)


# Calculate a pairwise association between all variables in a data-frame. In particular nominal vs nominal with Chi-square, numeric vs numeric with Pearson correlation, and nominal vs numeric with ANOVA.
# Adopted from https://stackoverflow.com/a/52557631/590437
mixed_assoc = function(df, cor_method="spearman", adjust_cramersv_bias=TRUE){
    df_comb = expand.grid(names(df), names(df),  stringsAsFactors = F) %>% set_names("X1", "X2")

    is_nominal = function(x) class(x) %in% c("factor", "character")
    # https://community.rstudio.com/t/why-is-purr-is-numeric-deprecated/3559
    # https://github.com/r-lib/rlang/issues/781
    is_numeric <- function(x) { is.integer(x) || is_double(x)}

    f = function(xName,yName) {
        x =  pull(df, xName)
        y =  pull(df, yName)

        result = if(is_nominal(x) && is_nominal(y)){
            # use bias corrected cramersV as described in https://rdrr.io/cran/rcompanion/man/cramerV.html
            cv = cramerV(as.character(x), as.character(y), bias.correct = adjust_cramersv_bias)
            data.frame(xName, yName, assoc=cv, type="cramersV")

        }else if(is_numeric(x) && is_numeric(y)){
            correlation = cor(x, y, method=cor_method, use="complete.obs")
            data.frame(xName, yName, assoc=correlation, type="correlation")

        }else if(is_numeric(x) && is_nominal(y)){
            # from https://stats.stackexchange.com/questions/119835/correlation-between-a-nominal-iv-and-a-continuous-dv-variable/124618#124618
            r_squared = summary(lm(x ~ y))$r.squared
            data.frame(xName, yName, assoc=sqrt(r_squared), type="anova")

        }else if(is_nominal(x) && is_numeric(y)){
            r_squared = summary(lm(y ~x))$r.squared
            data.frame(xName, yName, assoc=sqrt(r_squared), type="anova")

        }else {
            warning(paste("unmatched column type combination: ", class(x), class(y)))
        }

        # finally add complete obs number and ratio to table
        result %>% mutate(complete_obs_pairs=sum(!is.na(x) & !is.na(y)), complete_obs_ratio=complete_obs_pairs/length(x)) %>% rename(x=xName, y=yName)
    }

    # apply function to each variable combination
    map2_df(df_comb$X1, df_comb$X2, f)
}

Используя метод, мы можемлегко анализировать широкий спектр смешанных переменных фреймов данных:

mixed_assoc(iris)

              x            y      assoc        type complete_obs_pairs 
1  Sepal.Length Sepal.Length  1.0000000 correlation                150
2   Sepal.Width Sepal.Length -0.1667777 correlation                150
3  Petal.Length Sepal.Length  0.8818981 correlation                150
4   Petal.Width Sepal.Length  0.8342888 correlation                150
5       Species Sepal.Length  0.7865785       anova                150
6  Sepal.Length  Sepal.Width -0.1667777 correlation                150
7   Sepal.Width  Sepal.Width  1.0000000 correlation                150
25      Species      Species  1.0000000    cramersV                150

Это также можно использовать вместе с превосходным пакетом corrr, например, для построения корреляционной сетиграфик:

require(corrr)

msleep %>%
    select(- name) %>%
    mixed_assoc() %>%
    select(x, y, assoc) %>%
    spread(y, assoc) %>%
    column_to_rownames("x") %>%
    as.matrix %>%
    as_cordf %>%
    network_plot()

Answer 3

Вот решение tidyverse:

# example dataframe
df <- data.frame(
  group = c('A', 'A', 'A', 'A', 'A', 'B', 'C'),
  student = c('01', '01', '01', '02', '02', '01', '02'),
  exam_pass = c('Y', 'N', 'Y', 'N', 'Y', 'Y', 'N'),
  subject = c('Math', 'Science', 'Japanese', 'Math', 'Science', 'Japanese', 'Math')
) 

library(tidyverse)
library(lsr)

# function to get chi square p value and Cramers V
f = function(x,y) {
    tbl = df %>% select(x,y) %>% table()
    chisq_pval = round(chisq.test(tbl)$p.value, 4)
    cramV = round(cramersV(tbl), 4) 
    data.frame(x, y, chisq_pval, cramV) }

# create unique combinations of column names
# sorting will help getting a better plot (upper triangular)
df_comb = data.frame(t(combn(sort(names(df)), 2)), stringsAsFactors = F)

# apply function to each variable combination
df_res = map2_df(df_comb$X1, df_comb$X2, f)

# plot results
df_res %>%
  ggplot(aes(x,y,fill=chisq_pval))+
  geom_tile()+
  geom_text(aes(x,y,label=cramV))+
  scale_fill_gradient(low="red", high="yellow")+
  theme_classic()

Обратите внимание, что я использую пакет lsr для вычисления Cramers V с использованиемcramersV функция.