## -----------------------------------------------------------------------------
data("vert_df", package = "pciR")

head(vert_df)

## -----------------------------------------------------------------------------
vert_df$log_body_mass <- log(vert_df$body_mass)
vert_df$inv_range_area <- 1/vert_df$range_area
vert_df$inv_brood_size <- 1/vert_df$brood_size
vert_df$inv_protected_area <- 1/((vert_df$protected_area*vert_df$range_area+0.0001)) #add small number to avoid division by zero

## -----------------------------------------------------------------------------
traits_vertebrates <-
   vert_df[c("log_body_mass",
             "inv_range_area",
             "inv_brood_size",
             "inv_protected_area",
             "AHI")]

## -----------------------------------------------------------------------------
threats_2050_245 <-
   vert_df[c("clim_2050_245",
             "landuse_2050_245",
             "popdens_2050_245",
             "inv_threat")]

threats_2050_585 <-
   vert_df[c("clim_2050_585",
             "landuse_2050_585",
             "popdens_2050_585",
             "inv_threat")]

threats_2100_245 <-
   vert_df[c("clim_2100_245",
             "landuse_2100_245",
             "popdens_2100_245",
             "inv_threat")]

threats_2100_585 <-
   vert_df[c("clim_2100_585",
             "landuse_2100_585",
             "popdens_2100_585",
             "inv_threat")]

## -----------------------------------------------------------------------------
vert_pci_2050_45 <-
  pciR::pci(sp = vert_df$binomial,
           var_out = threats_2050_245,
           var_in = traits_vertebrates)

head(vert_pci_2050_45[order(vert_pci_2050_45$rank),])

## -----------------------------------------------------------------------------
vert_pci_2050_85 <-
  pciR::pci(sp = vert_df$binomial,
           var_out = threats_2050_585,
           var_in = traits_vertebrates)

head(vert_pci_2050_85[order(vert_pci_2050_85$rank),])

## -----------------------------------------------------------------------------
vert_pci_2100_45 <-
  pciR::pci(sp = vert_df$binomial,
           var_out = threats_2100_245,
           var_in = traits_vertebrates)

head(vert_pci_2100_45[order(vert_pci_2100_45$rank),])

## -----------------------------------------------------------------------------
vert_pci_2100_85 <-
  pciR::pci(sp = vert_df$binomial,
           var_out = threats_2100_585,
           var_in = traits_vertebrates)

head(vert_pci_2100_85[order(vert_pci_2100_85$rank),])

## -----------------------------------------------------------------------------
names(threats_2100_585)
# [1] "clim_2100_245"    "landuse_2100_245" "popdens_2100_245" "inv_threat"      

weights_threats <- c(3, 2, 1, 0)

## -----------------------------------------------------------------------------
names(traits_vertebrates)
# [1] "log_body_mass"      "inv_range_area"     "inv_brood_size"     "inv_protected_area" "AHI"      

weights_traits <-
matrix(1, 
       ncol(threats_2100_585), 
       ncol(traits_vertebrates))

rownames(weights_traits) <- names(threats_2100_585)
colnames(weights_traits) <- names(traits_vertebrates)

weights_traits

## -----------------------------------------------------------------------------
# climate change, all rows
weights_traits[,1] <- 2

# land use change, first two rows
weights_traits[1,2] <- 3
weights_traits[2,2] <- 3

weights_traits

## -----------------------------------------------------------------------------
vert_pci_2100_85_weighted <-
  pciR::pci(sp = vert_df$binomial,
           var_out = threats_2100_585,
           var_in = traits_vertebrates,
           weight_out = weights_threats,
           weight_in = weights_traits)

head(vert_pci_2100_85_weighted[order(vert_pci_2100_85_weighted$rank),])

## -----------------------------------------------------------------------------
rept_df <-
vert_df[vert_df$class == "Reptile",]

rept_df_iucn <-
  rept_df[rept_df$iucn_cat %in% 
            c("LC", "NT", "VU", "EN", "CR"),]

rept_df_iucn$iucn_cat_num <- 
  as.numeric(factor(rept_df_iucn$iucn_cat, 
         levels = c("LC", "NT", "VU", "EN", "CR")))

## -----------------------------------------------------------------------------
threats_2100_585_rept_iucn <- 
threats_2100_585[vert_df$binomial %in% rept_df_iucn$binomial,]

traits_rept_iucn <-
traits_vertebrates[vert_df$binomial %in% rept_df_iucn$binomial,]

## -----------------------------------------------------------------------------
optim_weights_rept <-
  pciR::optim_weights(sp = rept_df_iucn$binomial,
                     var_out = threats_2100_585_rept_iucn,
                     var_in = traits_rept_iucn,
                     reference = rept_df_iucn$iucn_cat_num,
                     type = "both",
                     control = list(maxit = 5))

## -----------------------------------------------------------------------------
optim_weights_rept$weight_out
optim_weights_rept$weight_in

## -----------------------------------------------------------------------------
threats_2100_585_rept <- 
threats_2100_585[vert_df$class == "Reptile",]

traits_rept <-
traits_vertebrates[vert_df$class == "Reptile",]

rept_pci_2100_85_opt <-
  pciR::pci(sp = rept_df$binomial,
           var_out = threats_2100_585_rept,
           var_in = traits_rept,
           weight_out = optim_weights_rept$weight_out,
           weight_in = optim_weights_rept$weight_in)

head(rept_pci_2100_85_opt[order(rept_pci_2100_85_opt$rank),])

