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

[Feature]: Do we need to put bias adjustment fraction within the likelihood? #185

Closed
1 task done
Bai-Li-NOAA opened this issue Jun 28, 2022 · 37 comments
Closed
1 task done
Labels
kind: question Further information is requested P1 high priority task
Milestone

Comments

@Bai-Li-NOAA
Copy link
Contributor

Bai-Li-NOAA commented Jun 28, 2022

Is your feature request related to a problem? Please describe.

While developing the recruitment likelihood component, the implementation subgroup came across a question that asked if we need the bias adjustment fraction within the likelihood. Can't we just use the standard likelihood formulation without modification?

The standard negative log-likelihood (nll) is

  1. $nll = \frac{n}{2}ln(2\pi)+\frac{1}{2}{\sum_{y=1}}^n(\frac{r_{y}^2}{\sigma_{R}^2}+ln(\sigma_{R}^2))$

where $n, r_{y}, and \sigma_{R}$ are the number of data points, recruitment deviations in log space, and standard deviation in log space respectively.

The nll that includes bias adjustment fraction based on A.3.10 from Methot and Wetzel 2013 is

  1. $nll = \frac{n}{2}ln(2\pi)+\frac{1}{2}{\sum_{y=1}}^n(\frac{r_{y}^2}{\sigma_{R}^2}+b_{y}ln(\sigma_{R}^2))$
    where $b_{y}$ is the bias adjustment fraction.

But we are not sure about back-transforming the equation 2 out to the original likelihood form and how the derivation justified:

  1. $\frac{1}{\sqrt{2\pi}^n {\sqrt{\sigma_{R}^2}}^{\sum{b_{y}}}}e^{-{\prod_{y=1}}^n \frac{r_{y}^2}{2\sigma_{R}^2}}$

Currently, we have implemented the standard nll (equation 1) in FIMS for M1. Do we need to modify the code and include bias adjustment fraction?

@KyleShertzer-NOAA, @iantaylor-NOAA, @Andrea-Havron-NOAA, please feel free to edit the description.

Describe the solution you would like.

NA

Describe alternatives you have considered

NA

Statistical validity, if applicable

No response

Describe if this is needed for a management application

No response

Additional context

No response

Code of Conduct

  • I agree to follow this project's Code of Conduct
@Bai-Li-NOAA Bai-Li-NOAA added the status: triage_needed this is not approved for this milestone, do not work on it yet label Jun 28, 2022
@Bai-Li-NOAA Bai-Li-NOAA added the kind: question Further information is requested label Jun 28, 2022
@iantaylor-NOAA
Copy link
Contributor

Explanation of bias adjustment in equation (2) from Methot and Taylor (2011) is as follows. This is less a formal statistical derivation and more a pragmatic solution that provides for desirable properties. @Rick-Methot-NOAA, please feel free to add more info on this.
image
image

@ChristineStawitz-NOAA
Copy link
Contributor

We are inclined to leave this as is and revisit post-M1 since existing models without the bias adjustment fraction in the likelihood (AMAK, MAS, BAM, ASAP) seem to yield comparable results.

@kellijohnson-NOAA mentions this likely only matters if you are missing many years of age comp data and have only a couple years in the beginning that aren't modeled fully and the rest of the years $b_y$ is close to 1

@ChristineStawitz-NOAA ChristineStawitz-NOAA removed the status: triage_needed this is not approved for this milestone, do not work on it yet label Jun 28, 2022
@ChristineStawitz-NOAA ChristineStawitz-NOAA added this to the MQ milestone Jun 28, 2022
@Rick-Methot-NOAA
Copy link

Agree with Kelli. When models are supplied with good age data for entire time series, the issue of bias adjustment is ignorable. Bias adjustment is most important when age data are weak or missing. SS3 models that use only length data and still try to estimate recdevs are the situation most in need of bias adjustment concept.

@Andrea-Havron-NOAA Andrea-Havron-NOAA added status: needs discussion Dialogue is needed before a decision can be made. P1 high priority task labels Jul 11, 2023
@Cole-Monnahan-NOAA
Copy link
Contributor

@Rick-Methot-NOAA Do we even need to do this manual bias correction to recruitment when TMB has a generic form built in already? This is a penalized ML thing to do and we're not longer bound by ADMB's limitations.

https://doi.org/10.1016/j.fishres.2015.11.016

In this paper, we have demonstrated that a new epsilon estimator improves upon the bias-correction algorithm previously introduced by Methot and Taylor (2011) for models that treat recruitment deviations as random effects.

@Rick-Methot-NOAA
Copy link

Good question Cole. I hope you can avoid the bias adjustment for the reasons demonstrated in that paper by Jim and Kasper. The paper by Ian and me showed that the bias adjustment was unnecessary when using ADMB's MCMC and SS3 turns off the recruitment bias adjustment when running in MCMC mode.
I do not know if TMB will adequately deal with the situation that requires the bias adjustment ramp in ADMB's MLE. So, I recommend that someone try TMB with a ramp in data quality regarding recruitment and see if TMB adequately adjusts for that.

@Cole-Monnahan-NOAA
Copy link
Contributor

@James-Thorson-NOAA could you chime in on this topic? The basic question is whether FIMS should build in the recruitment bias adjustment ramp or if your generic method totally supersedes that?

@Rick-Methot-NOAA The ramp option won't match the generic one, so how would we set up this experiment? Put the ramp into WHAM and compare the two? Compare it to MCMC?

@Andrea-Havron-NOAA
Copy link
Collaborator

@timjmiller is looking into bias adjustment as part of WHAM's research track and could potentially provide guidance on how we develop this in FIMS.

@Rick-Methot-NOAA
Copy link

I would:

  1. simulate a time series of recruitment at some sigmaR
  2. sample composition data from the generated population. Let the effN of that sampling increase over time
  3. run MLE and MCMC estimation to see if TMB's MLE can reproduce the arithmetic mean recruitment in the early years.

@Cole-Monnahan-NOAA
Copy link
Contributor

OK I'll try it with a simple WHAM model. I'm assuming something like effN=seq(5,100, len=nyrs) would suffice?

@timjmiller
Copy link
Contributor

timjmiller commented Aug 1, 2023 via email

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 1, 2023 via email

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 1, 2023 via email

@Rick-Methot-NOAA
Copy link

Tim - population abundance must be mean unbiased because the catch coming from that population is mean unbiased.

@Cole-Monnahan-NOAA
Copy link
Contributor

I worked up a simple simulation example using WHAM (script to reproduce but note you have to install a modified branch of WHAM b/c it doesn't ADREPORT recruits). The surprising result is that bias correction does almost nothing in this example, even if I turn the age comp data way down. sigmaR=1.35 so it's really big.

image

The key code modifications are

effN <- floor(seq(5,100, len=44))
effN[10:20] <- 1
input2$data$catch_Neff[,1] <- effN
input2$data$index_Neff <- cbind(effN,effN)*0+1

I left some early fishery age data so it could estimate initial age structure, but then have some years with only a single MN sample for the fishery. The two surveys have 1 in every year just so I don't have to change the code as much.

It's possible (likely?) I'm doing the simulation wrong and I'll check again tomorrow. @timjmiller can probably most efficiently check this if it's worth pursuing further.

@timjmiller
Copy link
Contributor

timjmiller commented Aug 2, 2023 via email

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 2, 2023 via email

@Cole-Monnahan-NOAA
Copy link
Contributor

@timjmiller my script certainly corroborates that finding. I did find a mistake I was taking mean of the posterior so I switched to median and now they all match very closely.

But I'm confused b/c I have some years with very little age data and sigmaR is big so the bias correction should be impactful in those cases. I also tried with bias_correct_pe=0 and it made no difference. So I'm at a bit of a loss there. @James-Thorson-NOAA is this unexpected from your perspective?

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 2, 2023 via email

@ChristineStawitz-NOAA
Copy link
Contributor

@James-Thorson-NOAA your image didn't come through - can you add it again please?

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 2, 2023 via email

@Rick-Methot-NOAA
Copy link

Rick-Methot-NOAA commented Aug 2, 2023

In Cole's example, the estimated recruitments are showing large changes between years, which means there is a lot more data informing them than many of the long time series situations that the bias ramp in SS3 is designed for.

A typical situation with sigmaR=0.4 (current Petrale sole assessment) has the se of each recent log recruitment down to 0.1, whereas the se of the early recruitments is equal to 0.4 because they are only informed by sigmaR. All the early recruitments have a MLE value of 0.0 and they show nil autocorrelation.
Info from Ian Taylor
image

@Cole-Monnahan-NOAA
Copy link
Contributor

I'm trying to produce a situation like Rick's using a packaged WHAM example but it is turning out to be difficult. I got something close with these estimated recruits:

image

but it is quite contrived. I also can't get MCMC to work b/c the early Fdevs wander off to +/- Inf for reasons I don't understand. But we can compare the MLE estimates with and without the bias correction, below are the first 15 years of estimates of recruits. So it's definitely doing more now, but I'm surprised to see that sometimes the correction decreases and sometimes increases the estimate of recruits.

Really we should find an SS3 model like Rick's example that has been ported to WHAM or some other TMB platform so that we can more directly test this effect. Does anyone have that lying around? I can do the legwork.

>    uncorrected corrected mle_se
1        26240     26240      0
2        15331     15104   6503
3         5252      5648   3680
4         9068      8252   5526
5        63903     66185 111040
6        12197      9759  13889
7        14055     13389  17324
8        17901     23222  25598
9        20123     16308  26755
10       20119     17816  27244
11       22049     27124  31902
12       24268     23625  34294
13       23718     18418  31945
14       18823     25157  25127
15       17124     14196  15839

@James-Thorson-NOAA
Copy link

My memory is that a negative correction is not possible using Methot-Taylor, but does in fact happen using epsilon correction (and make sense) when the covariance is non-zero (and I think highly negative) between rec-devs.

but I agree that it would be easier to find a model that already has a non-trivial estimate of the bias-ramp, if the goal is to investigate simple ways of replacing the bias-ramp using options in TMB.

@Rick-Methot-NOAA
Copy link

Rick-Methot-NOAA commented Aug 2, 2023

Sorry Cole, I cannot help with porting to WHAM.

While on this bias adjustment topic, a side issue is that large recruitments are estimated more precisely than weak recruitments (see trend in the green values). Also note that sometimes the se on the early recruitments is >sigmaR. This might be consequence of sum-to-zero constraint. The lower se on large recruitments might be due to: (a) use of multinomial logL; (b) data are length comps and age comps with ageing error, so small recruitments disappear into the noise. A corollary of this as noted in M&T(2011) is that large positive recdevs will get estimated, but -2sigmaR devs don't get estimated because they are lost in the noise and the sigmaR penalty draws them towards 0.0.

image

@timjmiller
Copy link
Contributor

timjmiller commented Aug 3, 2023 via email

@Rick-Methot-NOAA
Copy link

Tim - no trend in recruitment. And Ian & I saw the same pattern in our simulation work in 2011 which had a much larger sample size than one example time series.

@Cole-Monnahan-NOAA
Copy link
Contributor

I think we're at a bit of a standstill on this issue. To really investigate this we need a more thorough example that uses the ramp in SS3 and is ported to WHAM or another TMB platform. Maybe @kellijohnson-NOAA knows of one?

If not, I suggest that we not include the bias adjustment ramp in FIMS. Thorson and Kristensen provide a theoretical argument and examples that the generic bias adjustment is better. It's also not used in other platforms. People bridging from SS3 can rerun their models without the ramp to aid that process, then turn bias correction on once in FIMS. If it turns out the ramp is necessary then we can revisit it.

@Rick-Methot-NOAA
Copy link

Much prefer to see your TMB example completed first. Can you rerun it with zero data for the first 10 years of the time series? Are estimated recruitments during those years mean or median unbiased?

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 9, 2023 via email

@Rick-Methot-NOAA
Copy link

I do not think CCSRA will demonstrate what we need because it only has data in the final year. The ramp concept is only relevant in situations where we have a trend over time in data quality regarding some RE.

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 9, 2023 via email

@Rick-Methot-NOAA
Copy link

Thanks Jim. seems that suitable example can be created with CCSRA.

@Cole-Monnahan-NOAA
Copy link
Contributor

If we can agree that's a suitable example and it already shows improved performance compared to the ramp then I think we're safe to close this issue with the conclusion to not include any bias correction into FIMS and let that be done in TMB instead. Any remaining objections? If not I'll close this out by the end of the week.

@James-Thorson-NOAA
Copy link

James-Thorson-NOAA commented Aug 9, 2023 via email

@Rick-Methot-NOAA
Copy link

Looks good to proceed with the TMB approach.

@jimianelli
Copy link
Contributor

jimianelli commented Aug 9, 2023 via email

@Cole-Monnahan-NOAA
Copy link
Contributor

Closing this with conclusion to not include the bias adjustment in FIMS and let TMB do it externally. We can revisit once we have a handful of SS3 models ported to FIMS and can explore the impacts of the ramp more thoroughly.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
kind: question Further information is requested P1 high priority task
Projects
None yet
Development

No branches or pull requests

10 participants