knitr::opts_chunk$set(echo = TRUE)

Aligning eofs from multiple sources

so first get some patterns

n_modes <- 6

cera_patterns <- get_patterns(cera, k = n_modes, rotate = FALSE)
cesm_patterns <- get_patterns(filter(cesm, between(time, 1901, 2010)), k = n_modes, rotate = FALSE)
ccsm_patterns <- get_patterns(filter(ccsm, time >= 1901), k = n_modes, rotate = FALSE)
noaa_patterns <- get_patterns(filter(noaa, between(time, 1901, 2010)), k = n_modes, rotate = FALSE)

resample patterns to same resolution (probably the coarser of the two?)

cross correlate all the spatial patterns?

plot(cera_patterns$eofs)
plot(ccsm_patterns$eofs)
t1 <- as_tibble(cera_patterns$eofs) %>% 
  pivot_wider(names_from = PC, values_from = weight) %>% 
  select(-x,-y) %>% remove_missing()

t2 <- as_tibble(ccsm_patterns$eofs) %>% 
  pivot_wider(names_from = PC, values_from = weight) %>% 
  select(-x,-y) %>% remove_missing

cor(t1, t2, use = 'complete.obs')
f1 <- as_tibble(st_crop(get_correlation(cera, cera_patterns),states_wus) ) %>% 
  pivot_wider(names_from = PC, values_from = SWE) %>% 
  select(-x,-y) %>% remove_missing()

f2 <- as_tibble(st_crop(get_correlation(ccsm, ccsm_patterns),states_wus) ) %>% 
  pivot_wider(names_from = PC, values_from = SWE) %>% 
  select(-x,-y) %>% remove_missing

library(corrplot)
corrplot(cor(f1, f2, use = 'complete.obs'))
cor(f1, f2, use = 'complete.obs') %>% as_tibble(rownames = 'PC') %>% pivot_longer(-PC) %>%
  group_by(PC) %>%
  filter(abs(value) == max(abs(value)))

cor(f1, f2, use = 'complete.obs') %>% as_tibble(rownames = 'PC') %>% pivot_longer(-PC) %>%
  group_by(name) %>%
  filter(abs(value) == max(abs(value)))
cos.sim=function(ma, mb){
  mat=tcrossprod(ma, mb)
  t1=sqrt(apply(ma, 1, crossprod))
  t2=sqrt(apply(mb, 1, crossprod))
  mat / outer(t1,t2)
}
cos.sim(as.matrix(f1), as.matrix(f2))
proxy::simil(f1, f2, method = 'cosine', by_rows = FALSE) %>% corrplot()
proxy::simil(t1, t2, method = 'cosine', by_rows = FALSE) %>% corrplot()
psych::factor.congruence(f1,f2) %>% corrplot()
psych::factor.congruence(t1,t2) %>% corrplot()
# the above are identical! so the burt/tucker index of factor congruence is the cosine similarity
EFAtools::COMPARE(as.matrix(t1), as.matrix(t2), reorder = 'congruence', na.rm = TRUE)$diff_corres_cros
t1 <- as_tibble(cera_patterns$eofs) %>% 
  pivot_wider(names_from = PC, values_from = weight) %>% 
  select(-x,-y) %>% remove_missing()

t2 <- as_tibble(ccsm_patterns$eofs) %>% 
  pivot_wider(names_from = PC, values_from = weight) %>% 
  select(-x,-y) %>% remove_missing

psych::factor.congruence(t1,t2) %>%
  as_tibble(rownames = 'ref') %>%
  pivot_longer(-ref) %>%
  group_by(ref) %>%
  arrange(-abs(value), .by_group = TRUE) %>%
  mutate(fit = case_when(abs(value) >= 0.98 ~ 'Excellent',
                         abs(value) >= 0.92 ~ 'Good',
                         abs(value) >= 0.82 ~ 'Borderline',
                         abs(value) >= 0.68 ~ 'Poor',
                         TRUE ~ 'Terrible'))

  filter(abs(value) >= 0.68 | abs(value) == max(abs(value)))
pca <- prcomp(data_RSE[1:150,], scale. = TRUE)
pca$rotation[, 1:2, drop = FALSE]

pca$rotation[, 1:2, drop = FALSE] %>% # drop = FALSE preserves PC names when there's only 1 PC
      `%*%`(diag(pca$sdev, 2, 2))

pca$sdev^2
library(EFA.dimensions)
PCAoutput_1 <- PCA(data_RSE[1:150,],   Nfactors = 2,rotate = 'none', verbose=FALSE)

PCAoutput_2 <- PCA(data_RSE[151:300,], Nfactors = 2, rotate='PROMAX', verbose=FALSE)

t4 <- PROCRUSTES(target=as.matrix(t1), loadings=as.matrix(t2), 
           type = 'orthogonal', verbose=FALSE)
library(vegan)
procrustes(as.matrix(t1), as.matrix(t2), scale = FALSE) ->t3
t3$rotation;psych::factor.congruence(t1,t2) 
psych::factor.congruence(t1,t4$loadingsPROC) #%>% diag # these are the correct phi values post rotation!
cor(t1,t4$loadingsPROC) %>% corrplot()

so the problem with both correlation, pearson and spearman, and cosine similarities is that there is still mixing across the modes and the patterns aren't precisely matched! rotate signs as needed

flag patterns with poor matching

visualize everything



nick-gauthier/tidyEOF documentation built on July 21, 2023, 8:25 a.m.