Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

vignette for missing strata #17

Open
JaneSullivan-NOAA opened this issue Oct 5, 2023 · 2 comments
Open

vignette for missing strata #17

JaneSullivan-NOAA opened this issue Oct 5, 2023 · 2 comments

Comments

@JaneSullivan-NOAA
Copy link
Collaborator

Add a vignette on how to fill in missing strata using rema when developing a design-based index of abundance (e.g., for use in an age-structured model). Example: GOA dover sole

@Cole-Monnahan-NOAA, in this case would be it be appropriate to use the predicted value SD of the missing strata (converted to the natural scale) when summed with the other design-based variance estimates? Also pls feel free to add thoughts or other ideas for vignettes!

@Cole-Monnahan-NOAA
Copy link

Yes I prefer that approach of filling in missing strata estimates with the rema output. A dover sole example would be welcome for sure.

@careymcg
Copy link

careymcg commented Oct 6, 2023

For a future vignette (though variance calculations need to be double-checked) - this code takes rema output for predictions and data used and stitches together a survey biomass index (and variance) using design-based estimates for every strata where an estimate exists (every strata that was sampled) and filling in with rema estimates for the missing strata:

`# stitch design based strata with rema where there are design-based gaps:
surv_years<-unique(new_data$Year)
rema_by_strata_2023<-output$biomass_by_strata %>%
mutate(obs_var = (obs_cv*obs)^2,pred=exp(log_pred),pred_var = (exp(log_pred) * (sqrt(exp(sd_log_pred^2)-1)))^2) %>%
select(c(year,strata,pred,pred_var,obs,obs_var)) %>%
filter(year %in% surv_years) %>%
mutate(index = coalesce(obs,pred),index_var = coalesce(obs_var,pred_var)) %>%
group_by(year)

index_total<-rema_by_strata_2023 %>% group_by(year) %>% summarise(tot_index = sum(index),tot_index_var = sum(index_var)) %>%
mutate(tot_index_cv = sqrt(tot_index_var)/tot_index) %>%
mutate(tot_index_log_se = sqrt(log(tot_index_cv^2+1)))`

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants