2
\$\begingroup\$

Have a look at this dataset. The"a" and "b" are used to make possible differentiate when the same variable was measured. In this case X1a and X1b access the same variable, but "a" was (suppose..) in the last year and "b" was this year.

Original dataframe
I just want to correlate "a" and "b" and plot it. Simple like that and I really imagine the following code can be improved.

The final plot will be this one: Correlation plot The data is fake, but it's virtually equal as the original dataset I'm working on.

 all_items <- data.frame("1a" = sample(1:5), 
                            "2a" = sample(1:5),
                            "3a" = sample(1:5),
                            "1b" = sample(1:5),
                            "2b" = sample(1:5),
                            "3b" = rep(sample(1:5),10))


    #matrix with correlation
    all_correlation <- cor(all_items, method = "spearman") %>% 
      as.data.frame()

    #filter
    all_correlation <- all_correlation %>% select(-c(ends_with("a"))) #columns

    #create a colum with the now name
    all_correlation <- all_correlation %>% 
      mutate(item = row.names(.)) %>% select(item, everything())

    #supress some rows
    all_correlation <- all_correlation %>%  filter(!grepl("b", item))

    #filter(stringr::str_detect(row.names(.), "b"))
    #get only the diagonal
    all_correlation <- data.frame(item=1:3,Result=diag(as.matrix(all_correlation[, -1])))

    #P Value
    all_correlation_p_value <- Hmisc::rcorr(as.matrix(all_items))$P %>% as.data.frame()

    #filter
    all_correlation_p_value <- all_correlation_p_value %>% select(-c(ends_with("a")))
    all_correlation_p_value <- all_correlation_p_value %>% mutate(item = row.names(.)) %>% select(item, everything())
    all_correlation_p_value <- all_correlation_p_value %>%  filter(!grepl("b", item))
    all_correlation_p_value <- data.frame(item=1:3,P_Valor=diag(as.matrix(all_correlation_p_value[, -1])))

    #General table with the correlation  results
    all_correlation <- right_join(all_correlation,all_correlation_p_value, by = "item")


#Plot
ggplot(all_correlation, aes(x=item, y=Result)) + 
  geom_point(aes(color=Result)) +
  geom_line() +
  annotate("text", x = all_correlation$item, 
           y=all_correlation$Result, 
           label = paste("P-value =",round(all_correlation$P_Valor,3)), hjust = -0.1, colour = "red") +
  scale_x_continuous(breaks = seq(1,3,1)) 
\$\endgroup\$
1
\$\begingroup\$

Here are two possible rewrites. For the first one, notice how Hmisc::rcorr already computes correlations ($r) and p-values ($P) matrices so you can base all your work from it, provided you extract the rights rows (a_vars) and columns (b_vars) and only keep the diagonal values:

correlations <- rcorr(as.matrix(all_items), type = "spearman")
a_vars <- ends_with("a", vars = names(all_items))
b_vars <- ends_with("b", vars = names(all_items))
all_correlation <- data.frame(
  item    = seq_along(a_vars),
  Result  = diag(correlations$r[a_vars, b_vars]),
  P_Valor = diag(correlations$P[a_vars, b_vars])
)

What I don't like too much about this approach is that it only uses 3 of the 36 correlations that were computed by rcorr. In this second approach, I start by splitting the data into two tables (one for a and one for b) so as to compute only 3 correlations via Map. I also switched from rcorr to the base cor.test which has (IMO) a more intuitive behavior: when given two vectors as input, it computes one correlation, not four.

a <- all_items %>% select(ends_with("a"))
b <- all_items %>% select(ends_with("b"))
cor_test <- Map(cor.test, a, b)
all_correlation <- data.frame(
  item    = seq_along(a),
  Result  = sapply(cor_test, `[[`, "estimate"),
  P_Valor = sapply(cor_test, `[[`, "p.value")
)
\$\endgroup\$
  • \$\begingroup\$ Wow!! Both are amazing solutions! Thanks much for replying my question. I will go with the first one because my knowledge of "sapply" still very incipient and (because of that) it is not easy to understand what sapply does. Thank you! \$\endgroup\$ – Luis Nov 25 '18 at 18:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.