seacode / gmacs

A generic size-structured stock assessment model
https://seacode.github.io/gmacs
19 stars 14 forks source link

size transition matrix calculation not consistent with molt increment calculation #68

Open wStockhausen opened 9 years ago

wStockhausen commented 9 years ago

"Current" code (develop branch I just cloned) in calc_growth_increments(): ''' for( h = 1; h <= nsex; h++ ) { for( l = 1; l <= nclass; l++ ) { molt_increment(h)(l) = alpha(h) - beta(h) * mid_points(l); } } ''' Note that molt_increment[] is based on the size bin MID_POINTS.

relevant code in calc_growth_transition(): for( h = 1; h <= nsex; h++ ) { At.initialize(); sbi = size_breaks / gscale(h);

    for( l = 1; l <= nclass; l++ )
    {
        dMeanSizeAfterMolt = (size_breaks(l) + molt_increment(h)(l)) / gscale(h);
        psi.initialize();
        for( ll = l; ll <= nclass+1; ll++ )
        {
            if( ll <= nclass+1 )
            {
                psi(ll) = cumd_gamma(sbi(ll),dMeanSizeAfterMolt);
            }
        }
        At(l)(l,nclass)  = first_difference(psi(l,nclass+1));
        At(l)(l,nclass)  = At(l)(l,nclass) / sum(At(l));
    }
    growth_transition(h) = At;

Note that dMeanSizeAfterMolt mixes size_breaks[l](size bin cutpoints) and molt_increment[h,l], which is based on size bin midpoints. This does not seem seems consistent. I believe the correct approach is to calculate dMeanSizeAfterMolt using the mid_points, then integrate the resulting pdf from lower to upper cutpoint for each "final" size bin. This would be equivalent to assuming ALL individuals in a size bin had the size at the midpoint before growth occurred. One could be equally consistent by basing molt_increment[h,l] on the lower cutpoint and use the code above in calc_growth_transition(), but this would be treating all individuals in a size bin as having the size at the lower cutpoint (which is not consistent with weight-at-size, etc., calculations). An alternative would be to subdivide the initial size bins (just for this calculation) into smaller sub-bins (might be good where size bins are quite large) and integrate over the sub-bins as well.

smartell commented 9 years ago

Hi Buck. I think I see your point where you want to integrate over all possible pre-molt sizes within each size interval. However, Would this not be completely confounded with the variance (gscale). And is that not what the transition probability reflects; both the initial pre-molt variation in size within the size interval, and the post-molt size across all size-intervals >= initialal size interval.