Skip to content

Commit

Permalink
update manual code
Browse files Browse the repository at this point in the history
  • Loading branch information
caitobrien committed Feb 19, 2025
1 parent b6497c7 commit b9595df
Showing 1 changed file with 28 additions and 15 deletions.
43 changes: 28 additions & 15 deletions data-raw/HydroSurvDOYTEMP_2_extract_model_predictions.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -25,11 +25,14 @@ Current code below is setup to manually run each df of choice through to extract
```{r data_import}
#import saved subset file
df<-read_csv(here("data-raw/HydroSurvDOYTEMP_data_subset", "natural_chinook_data_subset.csv"))
df<-read_csv(here("data-raw/HydroSurvDOYTEMP_data_subset", "hatchery_spsu_chinook_data_subset.csv"))
df<-read_csv(here("data-raw/HydroSurvDOYTEMP_data_subset", "hatchery_chinook_data_subset.csv"))
df<-read_csv(here("data-raw/HydroSurvDOYTEMP_data_subset", "natural_steelhead_data_subset.csv"))
df<-read_csv(here("data-raw/HydroSurvDOYTEMP_data_subset", "hatchery_steelhead_data_subset.csv"))
df<-read_csv(file.choose())
species<- "steelhead"
rear_type<-"natural"
```

## aggregate subset data
Expand Down Expand Up @@ -71,13 +74,13 @@ df.agg <- df%>%
mod_natural_chinook_doy <- readRDS(here::here("data-raw/models", "mod_natural_chinook_doy.rds"))
mod_natural_chinook_temp<- readRDS(here::here("data-raw/models", "mod_natural_chinook_temp.rds"))
mod_natural_steelhead_doy<- readRDS(here::here("Examples/results", "mod_natural_steelhead_doy.rds"))
mod_natural_steelhead_temp<- readRDS(here::here("Examples/results", "mod_natural_steelhead_temp.rds"))
mod_natural_steelhead_doy<- readRDS(here::here("data-raw/models", "mod_natural_steelhead_doy.rds"))
mod_natural_steelhead_temp<- readRDS(here::here("data-raw/models", "mod_natural_steelhead_temp.rds"))
mod_hatchery_chinook_doy<- readRDS(here::here("Examples/results", "mod_hatchery_chinook_doy.rds"))
mod_hatchery_chinook_temp<- readRDS(here::here("Examples/results", "mod_hatchery_chinook_temp.rds"))
mod_hatchery_steelhead_doy<- readRDS(here::here("Examples/results", "mod_hatchery_steelhead_doy.rds"))
mod_hatchery_steelhead_temp<- readRDS(here::here("Examples/results", "mod_hatchery_steelhead_temp.rds"))
mod_hatchery_chinook_doy<- readRDS(here::here("data-raw/models", "mod_hatchery_chinook_doy.rds"))
mod_hatchery_chinook_temp<- readRDS(here::here("data-raw/models", "mod_hatchery_chinook_temp.rds"))
mod_hatchery_steelhead_doy<- readRDS(here::here("data-raw/models", "mod_hatchery_steelhead_doy.rds"))
mod_hatchery_steelhead_temp<- readRDS(here::here("data-raw/models", "mod_hatchery_steelhead_temp.rds"))
```

Expand All @@ -91,18 +94,25 @@ scale.temp <- scale(as.numeric(df$BON.temp))
scale.x.temp <- attr(scale.temp, "scaled:center")
scale.sd.temp <- attr(scale.temp, "scaled:scale")
years <- if (species == "chinook") {
c(1993:1996, 1998:2018)
} else if (species == "steelhead" & rear_type == "natural") {
c(1994:2018) # Adjust this range if needed for natural steelhead
} else { # Assuming rear_type is "hatchery" for steelhead
c(1993:2018) # Adjust this range if needed for hatchery steelhead
}
newdata.doy <- expand_grid(
doyz = ((90:160) - scale.x.doy) / scale.sd.doy,
transport = c(0, 1),
year = c(1993:2018),
year = years,#c(1993:2018),
n = 1
)
newdata.temp <- expand_grid(
mean.tempz = seq_range(df$BON.tempz, n = 71), #to match DOY
transport = c(0, 1),
year = c(1993:2018),
year = years, #c(1993:2018),
n = 1
)
Expand All @@ -111,9 +121,9 @@ newdata<-bind_cols(newdata.doy,newdata.temp[,1])

```{r add_predictions_SAR_TI}
df.pred<-newdata %>%
tidybayes::add_linpred_draws(mod_natural_chinook_doy, re_formula = NULL, allow_new_levels = TRUE, ndraws = 10, transform =TRUE) %>%
tidybayes::add_epred_draws(mod_natural_steelhead_temp, re_formula = NULL, allow_new_levels = TRUE,transform =TRUE) %>%
tidybayes::median_qi() %>%
rename(SAR = .linpred,
rename(SAR = .epred,
SAR.lo = .lower,
SAR.hi = .upper) %>%
group_by(year, doyz) %>%
Expand All @@ -131,13 +141,13 @@ df.pred<-newdata %>%

```{r append_obsdata_to_predictions}
#repeat code for each to get prediction
df.pred.w.ch.doy<-df.pred %>%
df.pred.n.stl.temp<-df.pred %>%
mutate(year = as.character(year),
transport = as.character(transport)) %>%
left_join(select(df.agg, "transport", "year","doy", "sar.pit", "n"), by = c("transport", "year", "doy")) %>%
mutate( rear_type = "Natural-origin", #"Hatchery-origin",#
covariate = "Day-of-year (DOY)" , # "Temperature (°C)",#
species = "Chinook", #"Steelhead", #
mutate( rear_type = "Natural-origin",#"Hatchery-origin",#
covariate = "Temperature (°C)",#"Day-of-year (DOY)" ,#
species = "Steelhead", # "Chinook", #
transport = as.factor(transport),
year = as.factor(year)) %>%
rename( n.obs = n.y)
Expand All @@ -156,6 +166,9 @@ all<-rbind(df.pred.w.ch.doy, df.pred.w.ch.temp,
#saved file will be the dataset used within data for all plots
write.csv(all, row.names = F, here("data", "df_mod_predict.csv"))
df_mod_predict<-all
usethis::use_data(df_mod_predict, overwrite = TRUE)
```

Expand Down

0 comments on commit b9595df

Please sign in to comment.