Models of performance

knitr::opts_chunk$set(cache = TRUE, dev = c("png", "pdf"))
library(here)
source(here("src", "model_linear_grampians_validate_response.R"))

perf_southeast <- merge(
  subset(
    perf, ibra_subregion != "Greater Grampians" & !random_effect, -random_effect
  ),
  with(
    mds,
    setNames(
      aggregate(y, list(taxon, ibra_subregion), mean),
      c("taxon", "ibra_subregion", "prev")
    )
  )
)

perf_southeast$dist <- dists[perf_southeast$ibra_subregion]
perf_southeast$dist_env <- kldists[perf_southeast$ibra_subregion]
perf_southeast$dist_com <- community_dist[perf_southeast$ibra_subregion]

forms <- outer(
  c("log(AUC)", "log(`AUC-PR`)", "log(`Rel-AUC-PR`)", "deviance_explained"),
  c("dist", "dist_com", "dist_env"),
  function(perf, dist) {
    sprintf(
      "%s ~
        scale(%s) + scale(log10(prev)) + (1 | ibra_subregion) + (1 | taxon)",
      perf, dist
    )
  }
)
forms <- lapply(t(forms), as.formula)
models <- lapply(forms, lmer, data = perf_southeast)
model_summaries <- lapply(models, summary)

for(i in model_summaries) print(i)

Linear mixed model fit by REML ['lmerMod']
Formula: log(AUC) ~ scale(dist) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: -332

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.3428 -0.7379  0.0096  0.7661  1.9102 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.001877 0.04332 
 ibra_subregion (Intercept) 0.001759 0.04194 
 Residual                   0.028791 0.16968 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)        -0.410615   0.013550  -30.30
scale(dist)         0.016063   0.012649    1.27
scale(log10(prev)) -0.066578   0.008012   -8.31

Correlation of Fixed Effects:
            (Intr) scl(d)
scale(dist) -0.018       
scl(lg10())  0.011 -0.156
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(AUC) ~ scale(dist_com) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: -330.5

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.3552 -0.7389  0.0163  0.7768  1.9651 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.001837 0.04286 
 ibra_subregion (Intercept) 0.002051 0.04529 
 Residual                   0.028810 0.16973 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)        -0.410510   0.014140 -29.033
scale(dist_com)     0.003784   0.012535   0.302
scale(log10(prev)) -0.065402   0.007956  -8.221

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_cm) -0.045       
scl(lg10())  0.010 -0.078
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(AUC) ~ scale(dist_env) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: -330.4

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.3596 -0.7430  0.0210  0.7793  1.9767 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.001832 0.04281 
 ibra_subregion (Intercept) 0.002057 0.04535 
 Residual                   0.028814 0.16975 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)        -0.410346   0.014136 -29.029
scale(dist_env)     0.002127   0.012577   0.169
scale(log10(prev)) -0.065108   0.007959  -8.180

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_nv) -0.011       
scl(lg10())  0.006  0.083
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(`AUC-PR`) ~ scale(dist) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: 1306.2

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8022 -0.5834 -0.1157  0.4133  5.1971 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.15016  0.3875  
 ibra_subregion (Intercept) 0.04496  0.2120  
 Residual                   0.48988  0.6999  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                   Estimate Std. Error t value
(Intercept)        -3.55332    0.07498 -47.390
scale(dist)         0.02089    0.06095   0.343
scale(log10(prev))  1.74817    0.03533  49.478

Correlation of Fixed Effects:
            (Intr) scl(d)
scale(dist) -0.022       
scl(lg10())  0.034 -0.135
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(`AUC-PR`) ~ scale(dist_com) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: 1306.4

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8098 -0.5849 -0.1167  0.4142  5.1958 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.15038  0.3878  
 ibra_subregion (Intercept) 0.04535  0.2130  
 Residual                   0.48983  0.6999  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                     Estimate Std. Error t value
(Intercept)        -3.5527588  0.0751905 -47.250
scale(dist_com)     0.0001471  0.0570058   0.003
scale(log10(prev))  1.7498628  0.0351037  49.848

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_cm) -0.041       
scl(lg10())  0.034 -0.071
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(`AUC-PR`) ~ scale(dist_env) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: 1304.1

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8500 -0.5779 -0.1164  0.3912  5.1967 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.15036  0.3878  
 ibra_subregion (Intercept) 0.03535  0.1880  
 Residual                   0.49020  0.7001  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                   Estimate Std. Error t value
(Intercept)        -3.55311    0.07126 -49.859
scale(dist_env)     0.08489    0.05242   1.619
scale(log10(prev))  1.75405    0.03503  50.077

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_nv) -0.004       
scl(lg10())  0.033  0.091
Linear mixed model fit by REML ['lmerMod']
Formula: 
log(`Rel-AUC-PR`) ~ scale(dist) + scale(log10(prev)) + (1 | ibra_subregion) +  
    (1 | taxon)
   Data: ..1

REML criterion at convergence: 1853.6

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-5.4799 -0.5305  0.0259  0.5563  3.5576 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.4852   0.6966  
 ibra_subregion (Intercept) 0.1177   0.3431  
 Residual                   1.2840   1.1332  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)        -4.212985   0.126246 -33.371
scale(dist)        -0.008537   0.098988  -0.086
scale(log10(prev))  1.831862   0.057676  31.762

Correlation of Fixed Effects:
            (Intr) scl(d)
scale(dist) -0.023       
scl(lg10())  0.036 -0.134
Linear mixed model fit by REML ['lmerMod']
Formula: log(`Rel-AUC-PR`) ~ scale(dist_com) + scale(log10(prev)) + (1 |  
    ibra_subregion) + (1 | taxon)
   Data: ..1

REML criterion at convergence: 1853.7

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-5.4815 -0.5301  0.0238  0.5570  3.5567 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.4848   0.6963  
 ibra_subregion (Intercept) 0.1173   0.3425  
 Residual                   1.2842   1.1332  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)        -4.212775   0.126211 -33.379
scale(dist_com)    -0.008373   0.092123  -0.091
scale(log10(prev))  1.831494   0.057300  31.963

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_cm) -0.040       
scl(lg10())  0.036 -0.071
Linear mixed model fit by REML ['lmerMod']
Formula: log(`Rel-AUC-PR`) ~ scale(dist_env) + scale(log10(prev)) + (1 |  
    ibra_subregion) + (1 | taxon)
   Data: ..1

REML criterion at convergence: 1849.7

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-5.4540 -0.5303  0.0178  0.5535  3.5575 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.48480  0.6963  
 ibra_subregion (Intercept) 0.07647  0.2765  
 Residual                   1.28551  1.1338  
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                   Estimate Std. Error t value
(Intercept)        -4.21362    0.11657 -36.148
scale(dist_env)     0.17497    0.08007   2.185
scale(log10(prev))  1.83788    0.05704  32.223

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_nv) -0.002       
scl(lg10())  0.037  0.097
Linear mixed model fit by REML ['lmerMod']
Formula: deviance_explained ~ scale(dist) + scale(log10(prev)) + (1 |  
    ibra_subregion) + (1 | taxon)
   Data: ..1

REML criterion at convergence: 681.5

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-9.4952 -0.3083  0.0558  0.4052  5.5121 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.169691 0.41194 
 ibra_subregion (Intercept) 0.002602 0.05101 
 Residual                   0.141361 0.37598 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                   Estimate Std. Error t value
(Intercept)         0.06073    0.05130   1.184
scale(dist)         0.02368    0.02308   1.026
scale(log10(prev)) -0.49763    0.01946 -25.570

Correlation of Fixed Effects:
            (Intr) scl(d)
scale(dist) -0.024       
scl(lg10())  0.048 -0.180
Linear mixed model fit by REML ['lmerMod']
Formula: deviance_explained ~ scale(dist_com) + scale(log10(prev)) + (1 |  
    ibra_subregion) + (1 | taxon)
   Data: ..1

REML criterion at convergence: 682.1

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-9.5068 -0.3073  0.0618  0.4155  5.4828 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.16796  0.40982 
 ibra_subregion (Intercept) 0.00252  0.05019 
 Residual                   0.14176  0.37651 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                   Estimate Std. Error t value
(Intercept)         0.06125    0.05104   1.200
scale(dist_com)     0.01608    0.02131   0.754
scale(log10(prev)) -0.49533    0.01924 -25.750

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_cm) -0.019       
scl(lg10())  0.046 -0.092
Linear mixed model fit by REML ['lmerMod']
Formula: deviance_explained ~ scale(dist_env) + scale(log10(prev)) + (1 |  
    ibra_subregion) + (1 | taxon)
   Data: ..1

REML criterion at convergence: 682.7

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-9.5209 -0.3061  0.0622  0.4137  5.4822 

Random effects:
 Groups         Name        Variance Std.Dev.
 taxon          (Intercept) 0.166756 0.40836 
 ibra_subregion (Intercept) 0.002947 0.05429 
 Residual                   0.141787 0.37655 
Number of obs: 560, groups:  taxon, 82; ibra_subregion, 18

Fixed effects:
                    Estimate Std. Error t value
(Intercept)         0.062073   0.051143   1.214
scale(dist_env)    -0.005124   0.020925  -0.245
scale(log10(prev)) -0.495066   0.019348 -25.587

Correlation of Fixed Effects:
            (Intr) scl(_)
scl(dst_nv) 0.007        
scl(lg10()) 0.044  0.124 

var_comps <- lapply(model_summaries, getElement, "varcor")
var_comps <- data.frame(
  rep(c("AUROC", "AUPRC", "Relative AUPRC", "Deviance Explained"), each = 3),
  c("Geographic", "Community", "Environmental"),
  sqrt(vapply(var_comps, getElement, 0, "taxon")),
  sqrt(vapply(var_comps, getElement, 0, "ibra_subregion")),
  vapply(var_comps, attr, 0, "sc")
)
names(var_comps) <-
  c("Performance metric", "Distance metric", "Taxon", "Region", "Residual")

vtab_perf <- xtable(
  subset(var_comps, `Performance metric` != "Relative AUPRC"),
  caption = "Variance components (on standard deviation scale) for models of
             performance metrics vs. three distance to test region metrics.",
  label = "tab:varcovperf",
  align = "lll|lll", digits = 2
)

dir.create(
  here("notebooks/models_of_performance_files/figure-html"),
  showWarnings = FALSE, recursive = TRUE
)

print(
  vtab_perf,
  file = here(
    "notebooks/models_of_performance_files/figure-html/vtab_perf.tex"
  ),
  include.rownames = FALSE,
  table.placement = "htbp!",
  caption.placement = "top",
  hline.after = c(-1, 0, nrow(vtab_perf)),
  sanitize.text.function = identity
)

fixefs_perf <- lapply(model_summaries, getElement, "coefficients")
fixefs_perf <- lapply(fixefs_perf, function(x) x[, 1:2])
fixefs_perf <- do.call(rbind, fixefs_perf)
fixefs_perf <- cbind(
  rep(c("AUROC", "AUPRC", "Relative AUPRC", "Deviance Explained"), each = 9),
  rep(c("Geographic", "Community", "Environmental"), each = 3),
  c("$\\beta_0$", "Distance", "Prevalence"),
  data.frame(fixefs_perf, row.names = NULL)
)
names(fixefs_perf) <- c(
  "Performance metric", "Distance metric", "Coefficient", "$\\mu$", "$\\sigma$"
)

fixefs_perf <- xtable(
  subset(fixefs_perf, `Performance metric` != "Relative AUPRC"),
  caption = "Fixed effects for models of
             performance metrics vs. three distance to test region metrics.",
  label = "tab:fixefsperf",
  align = "llll|ll", digits = 2
)

print(
  fixefs_perf,
  file = here(
    "notebooks/models_of_performance_files/figure-html/fixefs_perf.tex"
  ),
  include.rownames = FALSE,
  table.placement = "htbp!",
  caption.placement = "top",
  hline.after = c(-1, 0, nrow(fixefs_perf)),
  sanitize.text.function = identity
)

gg_and_se_taxa <- subset(mds, y)

gg_and_se_taxa <- names(
  Filter(
    function(x) x > 4,
    table(unique(gg_and_se_taxa[, c("taxon", "ibra_subregion")])$taxon)
  )
)

gg_and_se_taxa <- intersect(gg_and_se_taxa, unique(mdg$taxon))

perf_southeast$gg_and_se_taxa <- perf_southeast$taxon %in% gg_and_se_taxa

ggplot(perf_southeast) +
aes(prev, `AUC`, color = gg_and_se_taxa) +
geom_point() +
scale_color_grey() +
scale_x_log10() +
xlab("Prevalence") +
ylab("AUC Receiver-Operator") +
theme_bw() +
theme(legend.position = "none")

ggplot(perf_southeast) +
aes(`Rand-AUC-PR`, `AUC-PR`) +
geom_point(color = "darkgrey") +
xlab("AUC Precision-Recall of random classifier") +
ylab("AUC Precision-Recall") +
theme_bw()

ggplot(perf_southeast) +
aes(`Rand-AUC-PR`, `Rel-AUC-PR`) +
geom_point(color = "darkgrey") +
xlab("AUC Precision-Recall of random classifier") +
ylab("Relative AUC Precision-Recall") +
theme_bw()

ggplot(perf_southeast) +
aes(prev, `AUC-PR`) +
geom_point(color = "darkgrey") +
scale_x_log10()+
scale_y_log10() +
geom_line(
  aes(x = prev, y = 1 + (1 - prev) * log(1 - prev) / prev), color = "black",
  linetype = "dotted"
) +
xlab("Prevalence") +
ylab("AUC Precision-Recall") +
theme_bw()

pvm_data <- pivot_wider(
  responses,
  names_from = Var,
  values_from = c(Estimate, Std..Error, Prediction, Predicted..Error)
)

pvm_data <- filter(pvm_data, Region != "Greater Grampians")

pvm_data$dist_geo <- dists[pvm_data$Region]
pvm_data$dist_env <- kldists[pvm_data$Region]
pvm_data$dist_com <- community_dist[pvm_data$Region]

pvm_data <- mutate(
  pvm_data,
  mlq_mc = abs(Prediction_mlq - Estimate_mlq),
  twi_mc = abs(Prediction_twi - Estimate_twi),
  r1k_mc = abs(Prediction_r1k - Estimate_r1k),
  tn_mc  = abs(Prediction_tn - Estimate_tn)
)

pvm_auroc <- lmer(
  qlogis(AUC) ~
    scale(log(dist_geo)) +
    scale(log(dist_com)) +
    scale(log(dist_env)) +
    scale(log(mlq_mc)) +
    scale(log(twi_mc)) +
    scale(log(r1k_mc)) +
    scale(log(tn_mc)) +
    (1 | Taxon) +
    (1 | Region),
  data = pvm_data
)

summary(pvm_auroc)

Linear mixed model fit by REML ['lmerMod']
Formula: 
qlogis(AUC) ~ scale(log(dist_geo)) + scale(log(dist_com)) + scale(log(dist_env)) +  
    scale(log(mlq_mc)) + scale(log(twi_mc)) + scale(log(r1k_mc)) +  
    scale(log(tn_mc)) + (1 | Taxon) + (1 | Region)
   Data: pvm_data

REML criterion at convergence: 1384.4

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.0592 -0.5922 -0.2064  0.3506  5.4023 

Random effects:
 Groups   Name        Variance Std.Dev.
 Taxon    (Intercept) 0.08163  0.2857  
 Region   (Intercept) 0.03505  0.1872  
 Residual             0.58551  0.7652  
Number of obs: 560, groups:  Taxon, 82; Region, 18

Fixed effects:
                      Estimate Std. Error t value
(Intercept)           0.877001   0.065799  13.328
scale(log(dist_geo))  0.009232   0.074625   0.124
scale(log(dist_com)) -0.048411   0.073272  -0.661
scale(log(dist_env))  0.056034   0.061162   0.916
scale(log(mlq_mc))   -0.007920   0.036112  -0.219
scale(log(twi_mc))   -0.091729   0.036888  -2.487
scale(log(r1k_mc))   -0.039376   0.035976  -1.095
scale(log(tn_mc))    -0.023449   0.035176  -0.667

Correlation of Fixed Effects:
               (Intr) scl(lg(dst_g)) scl(lg(dst_c)) scl(lg(dst_n)) scl(lg(m_))
scl(lg(dst_g))  0.011                                                         
scl(lg(dst_c)) -0.041 -0.569                                                  
scl(lg(dst_n))  0.016 -0.116         -0.301                                   
scl(lg(m_))     0.001  0.014         -0.033         -0.037                    
scl(lg(tw_))    0.003  0.055          0.006         -0.036         -0.027     
scl(lg(1_))     0.002 -0.072          0.030          0.004         -0.004     
scl(lg(tn_))    0.016 -0.029          0.017          0.003         -0.022     
               scl(lg(tw_)) s((1_)
scl(lg(dst_g))                    
scl(lg(dst_c))                    
scl(lg(dst_n))                    
scl(lg(m_))                       
scl(lg(tw_))                      
scl(lg(1_))    -0.037             
scl(lg(tn_))   -0.088       -0.066

pvm_auprc <- lmer(
  log(`AUC-PR` / prevalence) ~
    scale(log(dist_geo)) +
    scale(log(dist_com)) +
    scale(log(dist_env)) +
    scale(log(mlq_mc)) +
    scale(log(twi_mc)) +
    scale(log(r1k_mc)) +
    scale(log(tn_mc)) +
    (1 | Taxon) +
    (1 | Region),
  data = pvm_data
)

summary(pvm_auprc)

Linear mixed model fit by REML ['lmerMod']
Formula: 
log(`AUC-PR`/prevalence) ~ scale(log(dist_geo)) + scale(log(dist_com)) +  
    scale(log(dist_env)) + scale(log(mlq_mc)) + scale(log(twi_mc)) +  
    scale(log(r1k_mc)) + scale(log(tn_mc)) + (1 | Taxon) + (1 |      Region)
   Data: pvm_data

REML criterion at convergence: 1285.2

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.6282 -0.5536 -0.1355  0.3694  5.7482 

Random effects:
 Groups   Name        Variance Std.Dev.
 Taxon    (Intercept) 0.16370  0.4046  
 Region   (Intercept) 0.03251  0.1803  
 Residual             0.45241  0.6726  
Number of obs: 560, groups:  Taxon, 82; Region, 18

Fixed effects:
                     Estimate Std. Error t value
(Intercept)           0.42741    0.07080   6.037
scale(log(dist_geo)) -0.04222    0.07029  -0.601
scale(log(dist_com)) -0.07632    0.06838  -1.116
scale(log(dist_env))  0.17035    0.05718   2.979
scale(log(mlq_mc))   -0.05378    0.03280  -1.640
scale(log(twi_mc))   -0.18672    0.03425  -5.452
scale(log(r1k_mc))   -0.14870    0.03289  -4.522
scale(log(tn_mc))    -0.04913    0.03183  -1.544

Correlation of Fixed Effects:
               (Intr) scl(lg(dst_g)) scl(lg(dst_c)) scl(lg(dst_n)) scl(lg(m_))
scl(lg(dst_g))  0.006                                                         
scl(lg(dst_c)) -0.035 -0.564                                                  
scl(lg(dst_n))  0.016 -0.115         -0.302                                   
scl(lg(m_))     0.000  0.012         -0.032         -0.035                    
scl(lg(tw_))    0.001  0.064          0.008         -0.038         -0.027     
scl(lg(1_))    -0.001 -0.063          0.031          0.001          0.018     
scl(lg(tn_))    0.020 -0.028          0.013          0.007         -0.021     
               scl(lg(tw_)) s((1_)
scl(lg(dst_g))                    
scl(lg(dst_c))                    
scl(lg(dst_n))                    
scl(lg(m_))                       
scl(lg(tw_))                      
scl(lg(1_))    -0.020             
scl(lg(tn_))   -0.097       -0.074

dir.create(
  here("notebooks/models_of_performance_files/figure-html"),
  showWarnings = FALSE, recursive = TRUE
)

feff_tab <- matrix("", 10, 6)
feff_tab[1, c(3, 5)] <- c("AUROC", "AUPRC/Prevalence")
feff_tab[2, 1:6] <-
  c("Fixed Effect", "", "$\\mu$", "$\\sigma$","$\\mu$", "$\\sigma$")
feff_tab[c(3, 4, 7), 1] <- c("Intercept", "Distance", "Error")
feff_tab[4:10, 2] <-
  c("Geographic", "Community", "Environmental", "MLQ", "TWI", "R1K", "TN")
feff_tab[3:10, 3:4] <-
  format(round(summary(pvm_auroc)$coefficients[, 1:2], 2), TRUE, 2)
feff_tab[3:10, 5:6] <-
  format(round(summary(pvm_auprc)$coefficients[, 1:2], 2), TRUE, 2)

feff_tab <- xtable(
  feff_tab, caption = "Fixed effects", label = "tab:feff",
  align = "lll|llll", digits = 1
)

print(
  feff_tab,
  file = here(
    "notebooks/models_of_performance_files/figure-html/fefftab.tex"
  ),
  include.rownames = FALSE,
  table.placement = "htbp!",
  caption.placement = "top",
  hline.after = c(-1, 2, nrow(feff_tab)),
  include.colnames = FALSE,
  sanitize.text.function = identity
)

reff_tab <- matrix(NA_character_, 4, 3)
reff_tab[1, 1:3] <- c("Random Effect", "AUROC", "AUPRC/Prevalence")
reff_tab[2:4, 1] <- c("Taxon", "Region", "Residual")
reff_tab[2:3, 2] <- format(sqrt(unlist(summary(pvm_auroc)$varcor)), digits = 1)
reff_tab[4, 2]   <- format(attr(summary(pvm_auroc)$varcor, "sc"), digits = 1)
reff_tab[2:3, 3] <- format(sqrt(unlist(summary(pvm_auprc)$varcor)), digits = 1)
reff_tab[4, 3]   <- format(attr(summary(pvm_auprc)$varcor, "sc"), digits = 1)

reff_tab <- xtable(
  reff_tab, caption = "Random effects", label = "tab:reff",
  align = "ll|ll", digits = 1
)

print(
  reff_tab,
  file = here(
    "notebooks/models_of_performance_files/figure-html/refftab.tex"
  ),
  include.rownames = FALSE,
  table.placement = "htbp!",
  caption.placement = "top",
  hline.after = c(-1, 1, nrow(reff_tab)),
  include.colnames = FALSE,
  sanitize.text.function = identity
)

ranef_traits <- unique(mds[c("taxon", "sla", "sm", "mh")])
ranef_traits$ranef <- ranef(pvm_auprc)$Taxon[ranef_traits$taxon, ]
ranef_traits <- pivot_longer(ranef_traits, sla:mh)

ggplot(ranef_traits) +
aes(value, ranef) +
geom_point() +
lims(x = c(-4, 4), y = c(-1.5, 1.5)) +
facet_grid(
  ~name,
  labeller = labeller(
    name = as_labeller(c(mh = "Maximum height", sla = "SLA", sm = "Seed Mass"))
  )
) +
xlab("Standard Deviations") +
ylab("Random effect") +
geom_text(
  data = summarise(group_by(ranef_traits, name), r = cor(ranef, value)),
  mapping = aes(
    x = -Inf, y = -Inf,
    label = sprintf("italic('r') == %s", round(r, digits = 2))
  ),
  hjust   = -0.1,
  vjust   = -1,
  parse = TRUE
) +
theme_bw() +
theme(strip.background = element_blank())