diff --git a/README.md b/README.md index a9e55825..618bc7a4 100644 --- a/README.md +++ b/README.md @@ -279,3 +279,6 @@ _____________________________________________________________ Last changed on: Source code: https://github.com/smartell/iSCAM _____________________________________________________________ + +## Code Changes + - March 30, 2015. Added time varying catchbility options (random walk). diff --git a/examples/makeproject b/examples/makeproject index 80d27455..c1e1005b 100755 --- a/examples/makeproject +++ b/examples/makeproject @@ -28,10 +28,10 @@ gsed -i "0,/^CTL=/s/^CTL=/CTL=$1/" $1/DATA/Makefile cp ../scripts/Template.dat $1/DATA/$1.dat cp ../scripts/Template.ctl $1/DATA/$1.ctl -cp ../scripts/Template.pfc $1/DATA/$1.pfc -cp ../scripts/Template.mpc $1/DATA/$1.mpc +cp ../scripts/Template.pfc $1/DATA/$1.pfc +cp ../scripts/Template.mpc $1/DATA/$1.mpc cp ../scripts/Template.scn $1/DATA/$1.scn -cp ../scripts/buildRdata.R $1/R/buildRdata.R +cp ../scripts/buildRdata.R $1/R/buildRdata.R cp ../scripts/collectAll.R $1/R/collectAll.R cp ../scripts/saveMSEdataframe.R $1/R/saveMSEdataframe.R diff --git a/fba/HALIBUT/DATA/DMVL1888/Halibut2014.ctl b/fba/HALIBUT/DATA/DMVL1888/Halibut2014.ctl index 271153cc..94e2037c 100644 --- a/fba/HALIBUT/DATA/DMVL1888/Halibut2014.ctl +++ b/fba/HALIBUT/DATA/DMVL1888/Halibut2014.ctl @@ -125,7 +125,7 @@ 2.0 # 3 -std in observed catches in first phase. 1.0 # 4 -std in observed catches in last phase. 1 # 5 -Assume unfished in first year (0=FALSE, 1=TRUE) -0.00 # 6 -Minimum proportion to consider in age-proportions for dmvlogistic +1.00 # 6 -Maternal effects power parameter (1.0 = no maternal effects) 0.30 # 7 -Mean fishing mortality for regularizing the estimates of Ft 0.10 # 8 -std in mean fishing mortality in first phase 2.00 # 9 -std in mean fishing mortality in last phase diff --git a/fba/HALIBUT/DATA/DMVL1988/Halibut2014.ctl b/fba/HALIBUT/DATA/DMVL1988/Halibut2014.ctl index 94eaed55..296fd100 100644 --- a/fba/HALIBUT/DATA/DMVL1988/Halibut2014.ctl +++ b/fba/HALIBUT/DATA/DMVL1988/Halibut2014.ctl @@ -1,9 +1,8 @@ -## ------------------------------------------------------------------------- ## -## CONTROL FILE TEMPLATE ## -## ------------------------------------------------------------------------- ## -## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## CONTROL FILE TEMPLATE ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## CONTROLS FOR LEADING PARAMETERS ## ## Prior descriptions: ## ## -0 uniform (0,0) ## @@ -11,23 +10,23 @@ ## -2 lognormal (p1=log(mu),p2=sig) ## ## -3 beta (p1=alpha,p2=beta) ## ## -4 gamma (p1=alpha,p2=beta) ## -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## npar 7 -## ival lb ub phz prior p1 p2 #parameter ## - 7.0 0.0 10 4 0 0.0 12.0 #log_ro ## - 0.75 0.2 1.0 4 3 3.00 2.00 #steepness ## - -1.89712 -3.0 2.0 -4 1 -1.203 0.15 #log_m g&b ## - 6.5 0.0 10 1 0 0.0 10 #log_avgrec ## - 7.5 0.0 10 1 0 0.0 10 #log_recinit ## - 0.00 0.00 1.00 -3 3 3.00 5.00 #rho ## - 0.6 0.01 5.0 -3 4 1.01 1.01 #sigma_r ## - -## ------------------------------------------------------------------------- ## +## ival lb ub phz prior p1 p2 #parameter ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + 7.0 0.0 10 4 0 0.0 10.0 #log_ro ## + 0.75 0.2 1.0 4 3 1.01 1.01 #steepness ## + -1.89712 -3.0 2.0 2 1 -1.89712 0.15 #log_m g&b ## + 6.5 0.0 10 1 0 0.0 10 #log_avgrec ## + 7.5 0.0 10 1 0 0.0 10 #log_recinit ## + 0.00 0.00 1.00 -3 3 3.00 5.00 #rho ## + 0.6 0.01 5.0 -3 4 1.01 1.01 #sigma_r ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## -## ------------------------------------------------------------------------- ## -## CONTROL PARAMETERS FOR AGE/SIZE COMPOSITION DATA FOR na_gears ## -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## CONTROL PARAMETERS FOR AGE/SIZE COMPOSITION DATA FOR na_gears ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## Likelihood type for each gear: ## -1 : multivariate logistic (dmvlogistic) ## -2 : multinomial, sample size based on input data @@ -35,20 +34,95 @@ ## -4 : logistic_normal, AR1 ## -5 : logistic_normal, AR2 ## -6 : multinomial with estimated effective sample size. -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## Number of columns == na_gears. - 1 6 6 ## : Gear Index - 6 6 6 ## : Likelihood type - 0.000 0.000 0.000 ## : Minimum proportion for aggregation & compression - 0.0000 0.0000 0.0000 ## : Small constant to add to comps & renormalize - -1 -1 -1 ## : phase for log_age_tau2 estimation. - -2 -2 -2 ## : phase for phi1 estimation : bounded (-1,1) AR1 - -2 -2 -2 ## : phase for phi2 estimation : bounded (0,1) AR2 - 2 2 2 ## : phase for degrees of freedom for student T. + 1 1 6 6 6 6 ## : Gear Index + 2 2 2 2 2 2 ## : Likelihood type + 0.000 0.000 0.000 0.000 0.000 0.000 ## : Minimum proportion for aggregation & compression + 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 ## : Small constant to add to comps & renormalize + -1 -1 -1 -1 -1 -1 ## : phase for log_age_tau2 estimation. + -2 -2 -2 -2 -2 -2 ## : phase for phi1 estimation : bounded (-1,1) AR1 + -2 -2 -2 -2 -2 -2 ## : phase for phi2 estimation : bounded (0,1) AR2 + -2 -2 -2 -2 -2 -2 ## : phase for degrees of freedom for student T. -12345 ## : int check (-12345) -## ------------------------------------------------------------------------- ## - +## ———————————————————————————————————————————————————————————————————————————————————— ## +## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## SELECTIVITY CONTROLS ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## - Each gear must have at least one selectivity and retention curve. +## - Use a -ve phase with the gear index to mirror another gear. Note +## that if you mirror another gear, it must have the same sel type and +## age and year nodes so that the arrays are the same shape & block years. ## +## • Index = gear index for selectivity curve. +## • sel_type = type of selectivity function (see Legend). +## • sel_mu = mean age/length 50% selectivity. +## • sel_sd = std in 50% selectivity +## • sex_dep = 0 -> no; 1 -> offset for sex 2. +## • size_nodes = # of nodes for age/size cubic spline. +## • year_nodes = # of nodes for time varying bicubic spline. +## • phz_mirror = phase of estimation (-ve phase to mirror selextivity index) +## • lam1 = penalty weight for 2nd differences (w = 1/(2•sig^2)). +## • lam2 = penalty weight for dome-shaped selectivity. +## • lam3 = penalty weight for time-varying selectivity. +## • start_block = year index for first year of selectivity curve. +## ———————————————————————————————————————————————————————————————————————————————————— ## +## sel_nBlocks ret_nBlocks ## Gear Description. + 1 1 ## Commercial retained + 1 1 ## Commercial discards + 1 1 ## Bycatch in non-directed fisheries. + 1 1 ## Sport + 1 1 ## Personal use + 2 1 ## Setline survey. +## ———————————————————————————————————————————————————————————————————————————————————— ## +## Selectivity P(capture of all size/age) +## slx_dControls +## • index for sex (0=both, 1=female, 2=male) +## sel sel sel age year phz or start end ## +## Index type mu sd sex nodes nodes mirror lam1 lam2 lam3 | block block ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + 1 5 10 2.0 0 5 6 2 12.5 2.0 50.5 1888 2014 + 2 1 3.0 1.5 0 0 0 -2 0.0 0.0 0.0 1888 2014 + 3 1 4.0 2.5 0 0 0 -3 0.0 0.0 0.0 1888 2014 + 4 1 4.0 2.0 0 0 0 -4 0.0 0.0 0.0 1888 2014 + 5 5 10 2.0 0 5 6 -1 0.0 0.0 0.0 1888 2014 + 6 5 3.0 2.0 1 5 6 2 12.5 2.0 50.5 1980 2014 + 6 5 3.0 2.0 2 5 6 2 12.5 2.0 50.5 1980 2014 +## ———————————————————————————————————————————————————————————————————————————————————— ## +## Retention P(retaining size/age) +## ret_dControls +## • index for sex (0=both, 1=female, 2=male) +## sel sel sel age year phz or start end ## +## Index type mu sd sex nodes nodes mirror lam1 lam2 lam3 | block block ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + -1 3 10 2.0 1 5 0 -2 0.0 0.0 0.0 1888 2014 + -2 1 3.0 1.5 0 5 0 -2 0.0 2.0 50.5 1888 2014 + -3 1 4.0 2.5 0 5 0 -3 0.0 0.0 0.0 1888 2014 + -4 1 4.0 2.0 0 5 0 -4 0.0 0.0 0.0 1888 2014 + -5 3 10 2.0 1 5 0 -1 0.0 0.0 0.0 1980 2014 + -6 5 3.0 2.0 1 5 5 -2 12.5 2.0 50.5 1980 2014 + +## ———————————————————————————————————————————————————————————————————————————————————— ## +## LEGEND FOR SELECTIVITY TYPES (sel_type) ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## sel | No. | +## type | parameters | Description +## ———————————————————————————————————————————————————————————————————————————————————— ## +## 1 | 2 | • Logistic curve with mean and standard dev at age = p(50%) +## 2 | (nages-1) | • Age-specific selectivity coefficients for (sage):(nage-1) +## 3 | age_nodes | • Age-specific coefficients based on cubic-spline interpolation +## 4 | n*age_nodes| • Annual age-specific coeffs using cubic-spline interpolation +## 5 | nyr*nage | • Bicubic spline interpolation over time & age. +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + + + +## TO BE DEPRECATED ## ------------------------------------------------------------------------- ## ## SELECTIVITY PARAMETERS Columns for gear ## ## OPTIONS FOR SELECTIVITY (isel_type): ## @@ -65,11 +139,11 @@ ## sig=0.05 0.10 0.15 0.20 0.30 0.40 0.50 ## ## wt =200. 50.0 22.2 12.5 5.56 3.12 2.00 ## ## ------------------------------------------------------------------------- ## - 3 1 1 1 3 11 # 1 -selectivity type ivector(isel_type) for gear + 3 1 1 1 3 5 # 1 -selectivity type ivector(isel_type) for gear 82 3.0 4.0 4.0 3.0 65.0 # 2 -Age/length at 50% selectivity (logistic) 5.5 1.5 2.5 2.0 2.5 5.5 # 3 -STD at 50% selectivity (logistic) 5 0 0 0 5 5 # 4 -No. of age nodes for each gear (0=ignore) - 0 0 0 0 0 10 # 5 -No. of year nodes for 2d spline(0=ignore) + 0 0 0 0 0 5 # 5 -No. of year nodes for 2d spline(0=ignore) 2 -2 -3 -4 -1 2 # 6 -Phase of estimation (-1 for fixed) 0.0 0.0 0.0 0.0 0.0 12.5 # 7 -Penalty wt for 2nd differences w=1/(2*sig^2) 200.0 200.0 200.0 200.0 200.0 200.0 # 8 -Penalty wt for dome-shaped w=1/(2*sig^2) @@ -82,8 +156,8 @@ 1888 1888 1888 -1888 -## +1970 +## ## ———————————————————————————————————————————————————————————————————————————————————— ## ## TIME VARYING NATURAL MORTALIIY RATES ## ## ———————————————————————————————————————————————————————————————————————————————————— ## @@ -104,18 +178,20 @@ ## ———————————————————————————————————————————————————————————————————————————————————— ## ## ## -## ------------------------------------------------------------------------- ## -## PRIORS FOR SURVEY Q ## -## Prior type: ## -## 0 - uninformative prior ## -## 1 - normal prior density for log(q) ## -## 2 - random walk in q ## -## ------------------------------------------------------------------------- ## -2 # -number of surveys (nits) -0 0 # -prior type (see legend above -0 0 # -prior log(mean) -0 0 # -prior sd -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## ABUNDANCE OBSERVATION MODELS +## ———————————————————————————————————————————————————————————————————————————————————— ## +## QTYPE: +## 0 = FIXED SURVEY Q (specify log(mean) for prior log(mean)) +## 1 = CONSTANT Q (use MLE for q and optional informative prior) +## 2 = RANDOM WALK Q (use prior mean & sd for penalized random walk) +## ———————————————————————————————————————————————————————————————————————————————————— ## +2 # -number of surveys (n_it_nobs) + 2 2 # -QTYPE (see legend above) + 0 0 # -prior log(mean) + 0 0 # -prior sd (set to 0 for uniformative prior) + 1 1 # -Estimation phase +## ———————————————————————————————————————————————————————————————————————————————————— ## ## ## ------------------------------------------------------------------------- ## @@ -126,14 +202,14 @@ 1.0 # 3 -std in observed catches in first phase. 1.0 # 4 -std in observed catches in last phase. 0 # 5 -Assume unfished in first year (0=FALSE, 1=TRUE) -0.00 # 6 -Minimum proportion to consider in age-proportions for dmvlogistic -0.40 # 7 -Mean fishing mortality for regularizing the estimates of Ft +1.00 # 6 -Maternal effects power parameter (1.0 = no maternal effects) +0.30 # 7 -Mean fishing mortality for regularizing the estimates of Ft 0.10 # 8 -std in mean fishing mortality in first phase 2.00 # 9 -std in mean fishing mortality in last phase -3 # 10 -DEPRECATED phase for estimating m_deviations (use -1 to turn off mdevs) 0.1 # 11 -DEPRECATED std in deviations for natural mortality 12 # 12 -DEPRECATED number of estimated nodes for deviations in natural mortality -0.00 # 13 -fraction of total mortality that takes place prior to spawning +1.00 # 13 -fraction of total mortality that takes place prior to spawning 82 # 14 -number of perspective years to start assessment at. 0 # 15 -switch for IFD distribution in selectivity simulations ## diff --git a/fba/HALIBUT/DATA/Halibut2014.dat b/fba/HALIBUT/DATA/Halibut2014.dat index 42cd972d..967dfc0c 100644 --- a/fba/HALIBUT/DATA/Halibut2014.dat +++ b/fba/HALIBUT/DATA/Halibut2014.dat @@ -384,223 +384,223 @@ ## commercial_data: ## survey ## iyr index(it) gear area group sex log(SE) log(PE) timing - 1907 280 1 1 1 0 0.100 0.010 0.5 - 1910 271 1 1 1 0 0.100 0.010 0.5 - 1911 237 1 1 1 0 0.100 0.010 0.5 - 1912 176 1 1 1 0 0.100 0.010 0.5 - 1913 129 1 1 1 0 0.100 0.010 0.5 - 1914 124 1 1 1 0 0.100 0.010 0.5 - 1915 118 1 1 1 0 0.100 0.010 0.5 - 1916 137 1 1 1 0 0.100 0.010 0.5 - 1917 98 1 1 1 0 0.100 0.010 0.5 - 1918 96 1 1 1 0 0.100 0.010 0.5 - 1919 93 1 1 1 0 0.100 0.010 0.5 - 1920 96 1 1 1 0 0.100 0.010 0.5 - 1921 88 1 1 1 0 0.100 0.010 0.5 - 1922 73 1 1 1 0 0.100 0.010 0.5 - 1923 78 1 1 1 0 0.100 0.010 0.5 - 1924 74 1 1 1 0 0.100 0.010 0.5 - 1925 68 1 1 1 0 0.100 0.010 0.5 - 1926 67 1 1 1 0 0.100 0.010 0.5 - 1927 65 1 1 1 0 0.100 0.010 0.5 - 1928 58 1 1 1 0 0.100 0.010 0.5 - 1929 51 1 1 1 0 0.100 0.010 0.5 - 1930 46 1 1 1 0 0.100 0.010 0.5 - 1931 50 1 1 1 0 0.100 0.010 0.5 - 1932 60 1 1 1 0 0.100 0.010 0.5 - 1933 63 1 1 1 0 0.100 0.010 0.5 - 1934 62 1 1 1 0 0.100 0.010 0.5 - 1935 76 1 1 1 0 0.100 0.010 0.5 - 1936 71 1 1 1 0 0.100 0.010 0.5 - 1937 80 1 1 1 0 0.100 0.010 0.5 - 1938 88 1 1 1 0 0.100 0.010 0.5 - 1939 80 1 1 1 0 0.100 0.010 0.5 - 1940 81 1 1 1 0 0.100 0.010 0.5 - 1941 85 1 1 1 0 0.100 0.010 0.5 - 1942 90 1 1 1 0 0.100 0.010 0.5 - 1943 95 1 1 1 0 0.100 0.010 0.5 - 1944 110 1 1 1 0 0.100 0.010 0.5 - 1945 102 1 1 1 0 0.100 0.010 0.5 - 1946 101 1 1 1 0 0.100 0.010 0.5 - 1947 99 1 1 1 0 0.100 0.010 0.5 - 1948 99 1 1 1 0 0.100 0.010 0.5 - 1949 95 1 1 1 0 0.100 0.010 0.5 - 1950 95 1 1 1 0 0.100 0.010 0.5 - 1951 96 1 1 1 0 0.100 0.010 0.5 - 1952 110 1 1 1 0 0.100 0.010 0.5 - 1953 131 1 1 1 0 0.100 0.010 0.5 - 1954 133 1 1 1 0 0.100 0.010 0.5 - 1955 119 1 1 1 0 0.100 0.010 0.5 - 1956 129 1 1 1 0 0.100 0.010 0.5 - 1957 110 1 1 1 0 0.100 0.010 0.5 - 1958 121 1 1 1 0 0.100 0.010 0.5 - 1959 129 1 1 1 0 0.100 0.010 0.5 - 1960 132 1 1 1 0 0.100 0.010 0.5 - 1961 127 1 1 1 0 0.100 0.010 0.5 - 1962 115 1 1 1 0 0.100 0.010 0.5 - 1963 105 1 1 1 0 0.100 0.010 0.5 - 1964 100 1 1 1 0 0.100 0.010 0.5 - 1965 99 1 1 1 0 0.100 0.010 0.5 - 1966 100 1 1 1 0 0.100 0.010 0.5 - 1967 101 1 1 1 0 0.100 0.010 0.5 - 1968 103 1 1 1 0 0.100 0.010 0.5 - 1969 95 1 1 1 0 0.100 0.010 0.5 - 1970 91 1 1 1 0 0.100 0.010 0.5 - 1971 89 1 1 1 0 0.100 0.010 0.5 - 1972 78 1 1 1 0 0.100 0.010 0.5 - 1973 63 1 1 1 0 0.100 0.010 0.5 - 1974 61 1 1 1 0 0.100 0.010 0.5 - 1975 61 1 1 1 0 0.100 0.010 0.5 - 1976 55 1 1 1 0 0.100 0.010 0.5 - 1977 63 1 1 1 0 0.100 0.010 0.5 - 1978 71 1 1 1 0 0.100 0.010 0.5 - 1979 75 1 1 1 0 0.100 0.010 0.5 - 1980 94 1 1 1 0 0.100 0.010 0.5 - 1981 111 1 1 1 0 0.100 0.010 0.5 - 1982 127 1 1 1 0 0.100 0.010 0.5 - 1984 291 1 1 1 0 0.100 0.010 0.5 - 1985 357 1 1 1 0 0.036 0.010 0.5 - 1986 320 1 1 1 0 0.041 0.010 0.5 - 1987 321 1 1 1 0 0.040 0.010 0.5 - 1988 368 1 1 1 0 0.035 0.010 0.5 - 1989 358 1 1 1 0 0.025 0.010 0.5 - 1990 319 1 1 1 0 0.028 0.010 0.5 - 1991 318 1 1 1 0 0.038 0.010 0.5 - 1992 319 1 1 1 0 0.034 0.010 0.5 - 1993 373 1 1 1 0 0.099 0.010 0.5 - 1994 306 1 1 1 0 0.072 0.010 0.5 - 1995 330 1 1 1 0 0.036 0.010 0.5 - 1996 392 1 1 1 0 0.038 0.010 0.5 - 1997 405 1 1 1 0 0.025 0.010 0.5 - 1998 408 1 1 1 0 0.025 0.010 0.5 - 1999 393 1 1 1 0 0.023 0.010 0.5 - 2000 403 1 1 1 0 0.022 0.010 0.5 - 2001 362 1 1 1 0 0.041 0.010 0.5 - 2002 360 1 1 1 0 0.019 0.010 0.5 - 2003 329 1 1 1 0 0.018 0.010 0.5 - 2004 319 1 1 1 0 0.019 0.010 0.5 - 2005 296 1 1 1 0 0.017 0.010 0.5 - 2006 270 1 1 1 0 0.019 0.010 0.5 - 2007 252 1 1 1 0 0.020 0.010 0.5 - 2008 232 1 1 1 0 0.017 0.010 0.5 - 2009 222 1 1 1 0 0.018 0.010 0.5 - 2010 204 1 1 1 0 0.020 0.010 0.5 - 2011 198 1 1 1 0 0.015 0.010 0.5 - 2012 195 1 1 1 0 0.021 0.010 0.5 - 2013 179 1 1 1 0 0.017 0.010 0.5 - 2014 191 1 1 1 0 0.052 0.010 0.5 + 1907 280 1 1 1 0 0.100 0.050 0.5 + 1910 271 1 1 1 0 0.100 0.050 0.5 + 1911 237 1 1 1 0 0.100 0.050 0.5 + 1912 176 1 1 1 0 0.100 0.050 0.5 + 1913 129 1 1 1 0 0.100 0.050 0.5 + 1914 124 1 1 1 0 0.100 0.050 0.5 + 1915 118 1 1 1 0 0.100 0.050 0.5 + 1916 137 1 1 1 0 0.100 0.050 0.5 + 1917 98 1 1 1 0 0.100 0.050 0.5 + 1918 96 1 1 1 0 0.100 0.050 0.5 + 1919 93 1 1 1 0 0.100 0.050 0.5 + 1920 96 1 1 1 0 0.100 0.050 0.5 + 1921 88 1 1 1 0 0.100 0.050 0.5 + 1922 73 1 1 1 0 0.100 0.050 0.5 + 1923 78 1 1 1 0 0.100 0.050 0.5 + 1924 74 1 1 1 0 0.100 0.050 0.5 + 1925 68 1 1 1 0 0.100 0.050 0.5 + 1926 67 1 1 1 0 0.100 0.050 0.5 + 1927 65 1 1 1 0 0.100 0.050 0.5 + 1928 58 1 1 1 0 0.100 0.050 0.5 + 1929 51 1 1 1 0 0.100 0.050 0.5 + 1930 46 1 1 1 0 0.100 0.050 0.5 + 1931 50 1 1 1 0 0.100 0.050 0.5 + 1932 60 1 1 1 0 0.100 0.050 0.5 + 1933 63 1 1 1 0 0.100 0.050 0.5 + 1934 62 1 1 1 0 0.100 0.050 0.5 + 1935 76 1 1 1 0 0.100 0.050 0.5 + 1936 71 1 1 1 0 0.100 0.050 0.5 + 1937 80 1 1 1 0 0.100 0.050 0.5 + 1938 88 1 1 1 0 0.100 0.050 0.5 + 1939 80 1 1 1 0 0.100 0.050 0.5 + 1940 81 1 1 1 0 0.100 0.050 0.5 + 1941 85 1 1 1 0 0.100 0.050 0.5 + 1942 90 1 1 1 0 0.100 0.050 0.5 + 1943 95 1 1 1 0 0.100 0.050 0.5 + 1944 110 1 1 1 0 0.100 0.050 0.5 + 1945 102 1 1 1 0 0.100 0.050 0.5 + 1946 101 1 1 1 0 0.100 0.050 0.5 + 1947 99 1 1 1 0 0.100 0.050 0.5 + 1948 99 1 1 1 0 0.100 0.050 0.5 + 1949 95 1 1 1 0 0.100 0.050 0.5 + 1950 95 1 1 1 0 0.100 0.050 0.5 + 1951 96 1 1 1 0 0.100 0.050 0.5 + 1952 110 1 1 1 0 0.100 0.050 0.5 + 1953 131 1 1 1 0 0.100 0.050 0.5 + 1954 133 1 1 1 0 0.100 0.050 0.5 + 1955 119 1 1 1 0 0.100 0.050 0.5 + 1956 129 1 1 1 0 0.100 0.050 0.5 + 1957 110 1 1 1 0 0.100 0.050 0.5 + 1958 121 1 1 1 0 0.100 0.050 0.5 + 1959 129 1 1 1 0 0.100 0.050 0.5 + 1960 132 1 1 1 0 0.100 0.050 0.5 + 1961 127 1 1 1 0 0.100 0.050 0.5 + 1962 115 1 1 1 0 0.100 0.050 0.5 + 1963 105 1 1 1 0 0.100 0.050 0.5 + 1964 100 1 1 1 0 0.100 0.050 0.5 + 1965 99 1 1 1 0 0.100 0.050 0.5 + 1966 100 1 1 1 0 0.100 0.050 0.5 + 1967 101 1 1 1 0 0.100 0.050 0.5 + 1968 103 1 1 1 0 0.100 0.050 0.5 + 1969 95 1 1 1 0 0.100 0.050 0.5 + 1970 91 1 1 1 0 0.100 0.050 0.5 + 1971 89 1 1 1 0 0.100 0.050 0.5 + 1972 78 1 1 1 0 0.100 0.050 0.5 + 1973 63 1 1 1 0 0.100 0.050 0.5 + 1974 61 1 1 1 0 0.100 0.050 0.5 + 1975 61 1 1 1 0 0.100 0.050 0.5 + 1976 55 1 1 1 0 0.100 0.050 0.5 + 1977 63 1 1 1 0 0.100 0.050 0.5 + 1978 71 1 1 1 0 0.100 0.050 0.5 + 1979 75 1 1 1 0 0.100 0.050 0.5 + 1980 94 1 1 1 0 0.100 0.050 0.5 + 1981 111 1 1 1 0 0.100 0.050 0.5 + 1982 127 1 1 1 0 0.100 2.000 0.5 + 1984 291 1 1 1 0 0.100 0.050 0.5 + 1985 357 1 1 1 0 0.036 0.050 0.5 + 1986 320 1 1 1 0 0.041 0.050 0.5 + 1987 321 1 1 1 0 0.040 0.050 0.5 + 1988 368 1 1 1 0 0.035 0.050 0.5 + 1989 358 1 1 1 0 0.025 0.050 0.5 + 1990 319 1 1 1 0 0.028 0.050 0.5 + 1991 318 1 1 1 0 0.038 0.050 0.5 + 1992 319 1 1 1 0 0.034 0.050 0.5 + 1993 373 1 1 1 0 0.099 0.050 0.5 + 1994 306 1 1 1 0 0.072 0.050 0.5 + 1995 330 1 1 1 0 0.036 0.050 0.5 + 1996 392 1 1 1 0 0.038 0.050 0.5 + 1997 405 1 1 1 0 0.025 0.050 0.5 + 1998 408 1 1 1 0 0.025 0.050 0.5 + 1999 393 1 1 1 0 0.023 0.050 0.5 + 2000 403 1 1 1 0 0.022 0.050 0.5 + 2001 362 1 1 1 0 0.041 0.050 0.5 + 2002 360 1 1 1 0 0.019 0.050 0.5 + 2003 329 1 1 1 0 0.018 0.050 0.5 + 2004 319 1 1 1 0 0.019 0.050 0.5 + 2005 296 1 1 1 0 0.017 0.050 0.5 + 2006 270 1 1 1 0 0.019 0.050 0.5 + 2007 252 1 1 1 0 0.020 0.050 0.5 + 2008 232 1 1 1 0 0.017 0.050 0.5 + 2009 222 1 1 1 0 0.018 0.050 0.5 + 2010 204 1 1 1 0 0.020 0.050 0.5 + 2011 198 1 1 1 0 0.015 0.050 0.5 + 2012 195 1 1 1 0 0.021 0.050 0.5 + 2013 179 1 1 1 0 0.017 0.050 0.5 + 2014 191 1 1 1 0 0.052 0.050 0.5 ## survey_data: - 1977 1.4713 6 1 1 0 0.100 0.010 0.5 - 1978 1.1112 6 1 1 0 0.100 0.010 0.5 - 1980 2.0090 6 1 1 0 0.100 0.010 0.5 - 1981 2.6670 6 1 1 0 0.100 0.010 0.5 - 1982 2.8739 6 1 1 0 0.100 0.010 0.5 - 1983 2.8816 6 1 1 0 0.100 0.010 0.5 - 1984 9.3014 6 1 1 0 0.100 0.010 0.5 - 1985 8.9419 6 1 1 0 0.100 0.010 0.5 - 1986 6.2568 6 1 1 0 0.100 0.010 0.5 - 1996 12.8858 6 1 1 0 0.100 0.010 0.5 - 1997 7.7841 6 1 1 0 0.034 0.010 0.5 - 1998 6.9727 6 1 1 0 0.034 0.010 0.5 - 1999 6.2649 6 1 1 0 0.034 0.010 0.5 - 2000 6.6265 6 1 1 0 0.034 0.010 0.5 - 2001 6.3842 6 1 1 0 0.028 0.010 0.5 - 2002 6.3893 6 1 1 0 0.026 0.010 0.5 - 2003 5.7876 6 1 1 0 0.027 0.010 0.5 - 2004 6.4579 6 1 1 0 0.026 0.010 0.5 - 2005 5.9282 6 1 1 0 0.027 0.010 0.5 - 2006 5.4306 6 1 1 0 0.027 0.010 0.5 - 2007 5.7267 6 1 1 0 0.027 0.010 0.5 - 2008 5.6536 6 1 1 0 0.026 0.010 0.5 - 2009 5.4952 6 1 1 0 0.027 0.010 0.5 - 2010 5.1239 6 1 1 0 0.029 0.010 0.5 - 2011 5.0435 6 1 1 0 0.029 0.010 0.5 - 2012 5.5606 6 1 1 0 0.030 0.010 0.5 - 2013 4.7301 6 1 1 0 0.034 0.010 0.5 - 2014 5.1645 6 1 1 0 0.032 0.010 0.5 + 1977 1.4713 6 1 1 0 0.100 0.001 0.5 + 1978 1.1112 6 1 1 0 0.100 0.001 0.5 + 1980 2.0090 6 1 1 0 0.100 0.001 0.5 + 1981 2.6670 6 1 1 0 0.100 1.001 0.5 + 1982 2.8739 6 1 1 0 0.100 0.001 0.5 + 1983 2.8816 6 1 1 0 0.100 1.001 0.5 + 1984 9.3014 6 1 1 0 0.100 0.001 0.5 + 1985 8.9419 6 1 1 0 0.100 0.001 0.5 + 1986 6.2568 6 1 1 0 0.100 0.001 0.5 + 1996 12.8858 6 1 1 0 0.100 1.001 0.5 + 1997 7.7841 6 1 1 0 0.034 0.001 0.5 + 1998 6.9727 6 1 1 0 0.034 0.001 0.5 + 1999 6.2649 6 1 1 0 0.034 0.001 0.5 + 2000 6.6265 6 1 1 0 0.034 0.001 0.5 + 2001 6.3842 6 1 1 0 0.028 0.001 0.5 + 2002 6.3893 6 1 1 0 0.026 0.001 0.5 + 2003 5.7876 6 1 1 0 0.027 0.001 0.5 + 2004 6.4579 6 1 1 0 0.026 0.001 0.5 + 2005 5.9282 6 1 1 0 0.027 0.001 0.5 + 2006 5.4306 6 1 1 0 0.027 0.001 0.5 + 2007 5.7267 6 1 1 0 0.027 0.001 0.5 + 2008 5.6536 6 1 1 0 0.026 0.001 0.5 + 2009 5.4952 6 1 1 0 0.027 0.001 0.5 + 2010 5.1239 6 1 1 0 0.029 0.001 0.5 + 2011 5.0435 6 1 1 0 0.029 0.001 0.5 + 2012 5.5606 6 1 1 0 0.030 0.001 0.5 + 2013 4.7301 6 1 1 0 0.034 0.001 0.5 + 2014 5.1645 6 1 1 0 0.032 0.001 0.5 ## ## ------------------------------------------------------------------------- ## ## AGE COMPOSITION DATA (ROW YEAR, COL=AGE) Ragged object ## ## ------------------------------------------------------------------------- ## -3 # Number of gears with age-comps int(na_gears) -80 26 26 # Number of rows in the matrix ivector(na_gears) +6 # Number of gears with age-comps int(na_gears) +67 13 12 14 12 14 # Number of rows in the matrix ivector(na_gears) ## ivector(na_gears) of youngest age-class -6 6 6 +6 6 6 6 6 6 ## ivector(na_gears) of oldest age-class + group -25 25 25 +20 25 20 25 20 25 ## Effective sample size from iterative re-weighting (nscaler) -111.4450 28.7800 43.9589 +50 50 50 50 50 50 ## Age Flag (0=length comps, 1 = age comps) -1 1 1 +1 1 1 1 1 1 ## commercial catches ## year gear area group sex age_err | data columns (numbers or proportions) - 1935 1 1 1 0 4 0.1085 0.1645 0.1479 0.1896 0.1486 0.0762 0.0610 0.0375 0.0225 0.0165 0.0112 0.0055 0.0044 0.0021 0.0041 0.0000 0.0000 0.0000 0.0000 0.0000 - 1936 1 1 1 0 4 0.0729 0.1376 0.1219 0.1225 0.1496 0.1208 0.0864 0.0673 0.0388 0.0261 0.0183 0.0132 0.0099 0.0042 0.0105 0.0000 0.0000 0.0000 0.0000 0.0000 - 1937 1 1 1 0 4 0.0715 0.0929 0.1694 0.1244 0.1129 0.1208 0.1000 0.0802 0.0379 0.0364 0.0168 0.0108 0.0073 0.0060 0.0127 0.0000 0.0000 0.0000 0.0000 0.0000 - 1938 1 1 1 0 4 0.1129 0.1154 0.1649 0.1634 0.0741 0.0772 0.0835 0.0751 0.0559 0.0315 0.0184 0.0100 0.0049 0.0050 0.0080 0.0000 0.0000 0.0000 0.0000 0.0000 - 1939 1 1 1 0 4 0.0883 0.0863 0.1537 0.1421 0.1535 0.0809 0.0819 0.0690 0.0608 0.0323 0.0269 0.0116 0.0074 0.0010 0.0041 0.0000 0.0000 0.0000 0.0000 0.0000 - 1940 1 1 1 0 4 0.1354 0.0819 0.0997 0.1274 0.1618 0.1236 0.0779 0.0485 0.0484 0.0475 0.0202 0.0132 0.0081 0.0016 0.0049 0.0000 0.0000 0.0000 0.0000 0.0000 - 1941 1 1 1 0 4 0.1321 0.0840 0.1386 0.1137 0.1337 0.1125 0.1064 0.0596 0.0347 0.0352 0.0245 0.0127 0.0054 0.0022 0.0046 0.0000 0.0000 0.0000 0.0000 0.0000 - 1942 1 1 1 0 4 0.0924 0.1572 0.1021 0.1197 0.1161 0.1057 0.1015 0.0787 0.0431 0.0260 0.0257 0.0132 0.0100 0.0038 0.0047 0.0000 0.0000 0.0000 0.0000 0.0000 - 1943 1 1 1 0 4 0.0954 0.1839 0.2094 0.0969 0.0855 0.0744 0.0720 0.0754 0.0382 0.0258 0.0190 0.0104 0.0057 0.0036 0.0044 0.0000 0.0000 0.0000 0.0000 0.0000 - 1944 1 1 1 0 4 0.0616 0.0793 0.1749 0.1459 0.1004 0.0857 0.0946 0.0978 0.0637 0.0371 0.0230 0.0160 0.0071 0.0078 0.0051 0.0000 0.0000 0.0000 0.0000 0.0000 - 1945 1 1 1 0 4 0.0892 0.1572 0.1585 0.1452 0.0999 0.0767 0.0744 0.0780 0.0465 0.0273 0.0177 0.0127 0.0056 0.0055 0.0057 0.0000 0.0000 0.0000 0.0000 0.0000 - 1946 1 1 1 0 4 0.0457 0.1226 0.1913 0.1477 0.1363 0.0852 0.0767 0.0761 0.0457 0.0264 0.0177 0.0122 0.0056 0.0051 0.0057 0.0000 0.0000 0.0000 0.0000 0.0000 - 1947 1 1 1 0 4 0.0630 0.0716 0.1631 0.1503 0.1277 0.1112 0.0933 0.0852 0.0513 0.0297 0.0202 0.0150 0.0066 0.0058 0.0061 0.0000 0.0000 0.0000 0.0000 0.0000 - 1948 1 1 1 0 4 0.0630 0.0871 0.1647 0.1762 0.1370 0.0976 0.0798 0.0761 0.0457 0.0262 0.0181 0.0123 0.0055 0.0052 0.0054 0.0000 0.0000 0.0000 0.0000 0.0000 - 1949 1 1 1 0 4 0.0371 0.0577 0.1235 0.1477 0.1650 0.1137 0.1188 0.0934 0.0551 0.0264 0.0190 0.0173 0.0092 0.0058 0.0101 0.0000 0.0000 0.0000 0.0000 0.0000 - 1950 1 1 1 0 4 0.0190 0.0408 0.1129 0.1959 0.1516 0.1165 0.1108 0.0877 0.0771 0.0301 0.0239 0.0126 0.0063 0.0046 0.0101 0.0000 0.0000 0.0000 0.0000 0.0000 - 1951 1 1 1 0 4 0.0085 0.0519 0.1064 0.1568 0.1718 0.1291 0.1021 0.0955 0.0615 0.0518 0.0262 0.0169 0.0101 0.0041 0.0073 0.0000 0.0000 0.0000 0.0000 0.0000 - 1952 1 1 1 0 4 0.0226 0.0216 0.1210 0.1495 0.1397 0.1443 0.1017 0.0956 0.0746 0.0522 0.0331 0.0204 0.0094 0.0064 0.0079 0.0000 0.0000 0.0000 0.0000 0.0000 - 1953 1 1 1 0 4 0.0289 0.0584 0.0866 0.1703 0.1513 0.1355 0.1093 0.0797 0.0571 0.0451 0.0322 0.0181 0.0141 0.0039 0.0095 0.0000 0.0000 0.0000 0.0000 0.0000 - 1954 1 1 1 0 4 0.0387 0.0467 0.1034 0.1053 0.1854 0.1143 0.0980 0.0848 0.0646 0.0505 0.0429 0.0273 0.0173 0.0093 0.0116 0.0000 0.0000 0.0000 0.0000 0.0000 - 1955 1 1 1 0 4 0.0889 0.0744 0.0985 0.1275 0.1059 0.1523 0.1011 0.0768 0.0511 0.0404 0.0303 0.0231 0.0125 0.0079 0.0095 0.0000 0.0000 0.0000 0.0000 0.0000 - 1956 1 1 1 0 4 0.0940 0.1027 0.1119 0.0929 0.1271 0.0856 0.1181 0.0770 0.0622 0.0410 0.0257 0.0208 0.0151 0.0121 0.0139 0.0000 0.0000 0.0000 0.0000 0.0000 - 1957 1 1 1 0 4 0.0812 0.1069 0.1202 0.1188 0.0805 0.1076 0.0752 0.1030 0.0603 0.0487 0.0323 0.0228 0.0172 0.0097 0.0155 0.0000 0.0000 0.0000 0.0000 0.0000 - 1958 1 1 1 0 4 0.1450 0.1145 0.1120 0.1036 0.0985 0.0761 0.0772 0.0571 0.0762 0.0471 0.0288 0.0209 0.0163 0.0111 0.0154 0.0000 0.0000 0.0000 0.0000 0.0000 - 1959 1 1 1 0 4 0.1324 0.0669 0.1993 0.1381 0.1077 0.0862 0.0522 0.0593 0.0388 0.0387 0.0234 0.0180 0.0115 0.0105 0.0168 0.0000 0.0000 0.0000 0.0000 0.0000 - 1960 1 1 1 0 4 0.1126 0.0778 0.1170 0.2314 0.1287 0.0831 0.0570 0.0437 0.0374 0.0305 0.0256 0.0173 0.0120 0.0086 0.0171 0.0000 0.0000 0.0000 0.0000 0.0000 - 1961 1 1 1 0 4 0.0644 0.1012 0.1132 0.1301 0.2120 0.1110 0.0719 0.0492 0.0351 0.0305 0.0239 0.0200 0.0124 0.0087 0.0164 0.0000 0.0000 0.0000 0.0000 0.0000 - 1962 1 1 1 0 4 0.0655 0.0763 0.1229 0.1351 0.1461 0.1892 0.0738 0.0491 0.0333 0.0261 0.0252 0.0192 0.0131 0.0080 0.0172 0.0000 0.0000 0.0000 0.0000 0.0000 - 1963 1 1 1 0 4 0.0818 0.0735 0.1434 0.1511 0.1308 0.1204 0.1309 0.0533 0.0279 0.0179 0.0172 0.0150 0.0121 0.0086 0.0160 0.0000 0.0000 0.0000 0.0000 0.0000 - 1964 1 1 1 0 4 0.0504 0.0929 0.1002 0.1885 0.1384 0.1168 0.0987 0.0976 0.0381 0.0189 0.0146 0.0121 0.0086 0.0080 0.0161 0.0000 0.0000 0.0000 0.0000 0.0000 - 1965 1 1 1 0 4 0.0524 0.0817 0.1552 0.1160 0.1898 0.1161 0.0880 0.0725 0.0562 0.0209 0.0129 0.0093 0.0083 0.0058 0.0148 0.0000 0.0000 0.0000 0.0000 0.0000 - 1966 1 1 1 0 4 0.0471 0.0494 0.1096 0.1534 0.1143 0.1619 0.1010 0.0805 0.0661 0.0464 0.0217 0.0137 0.0106 0.0056 0.0187 0.0000 0.0000 0.0000 0.0000 0.0000 - 1967 1 1 1 0 4 0.0778 0.0711 0.0881 0.1643 0.1548 0.0963 0.1196 0.0687 0.0479 0.0374 0.0269 0.0136 0.0094 0.0067 0.0173 0.0000 0.0000 0.0000 0.0000 0.0000 - 1968 1 1 1 0 4 0.0777 0.1108 0.0990 0.1021 0.1611 0.1256 0.0794 0.0932 0.0489 0.0354 0.0268 0.0173 0.0080 0.0044 0.0102 0.0000 0.0000 0.0000 0.0000 0.0000 - 1969 1 1 1 0 4 0.1105 0.0591 0.1727 0.0987 0.0872 0.1288 0.1020 0.0659 0.0637 0.0329 0.0272 0.0203 0.0117 0.0066 0.0128 0.0000 0.0000 0.0000 0.0000 0.0000 - 1970 1 1 1 0 4 0.1088 0.0799 0.0987 0.1886 0.0915 0.0820 0.1064 0.0727 0.0506 0.0450 0.0242 0.0158 0.0122 0.0089 0.0147 0.0000 0.0000 0.0000 0.0000 0.0000 - 1971 1 1 1 0 4 0.1255 0.1017 0.1335 0.1164 0.1729 0.0706 0.0621 0.0688 0.0479 0.0387 0.0238 0.0147 0.0094 0.0061 0.0081 0.0000 0.0000 0.0000 0.0000 0.0000 - 1972 1 1 1 0 4 0.1139 0.1099 0.1348 0.1347 0.1125 0.1493 0.0582 0.0573 0.0440 0.0318 0.0209 0.0132 0.0065 0.0049 0.0082 0.0000 0.0000 0.0000 0.0000 0.0000 - 1973 1 1 1 0 4 0.0269 0.0661 0.1138 0.1376 0.1489 0.1187 0.1492 0.0621 0.0486 0.0480 0.0308 0.0160 0.0136 0.0054 0.0144 0.0000 0.0000 0.0000 0.0000 0.0000 - 1974 1 1 1 0 4 0.0202 0.0454 0.1028 0.1246 0.1325 0.1304 0.1035 0.1276 0.0472 0.0465 0.0485 0.0268 0.0151 0.0115 0.0173 0.0000 0.0000 0.0000 0.0000 0.0000 - 1975 1 1 1 0 4 0.0178 0.0423 0.1025 0.1322 0.1423 0.1332 0.1220 0.0778 0.0906 0.0343 0.0321 0.0265 0.0170 0.0080 0.0215 0.0000 0.0000 0.0000 0.0000 0.0000 - 1976 1 1 1 0 4 0.0491 0.0539 0.1244 0.1442 0.1340 0.1095 0.1065 0.0744 0.0602 0.0617 0.0210 0.0196 0.0148 0.0081 0.0185 0.0000 0.0000 0.0000 0.0000 0.0000 - 1977 1 1 1 0 4 0.0288 0.0583 0.1065 0.1362 0.1577 0.1079 0.1049 0.0848 0.0654 0.0491 0.0453 0.0172 0.0137 0.0075 0.0165 0.0000 0.0000 0.0000 0.0000 0.0000 - 1978 1 1 1 0 4 0.0390 0.0691 0.1275 0.1387 0.1465 0.1236 0.1018 0.0691 0.0590 0.0425 0.0320 0.0190 0.0103 0.0065 0.0153 0.0000 0.0000 0.0000 0.0000 0.0000 - 1979 1 1 1 0 4 0.0165 0.0644 0.1065 0.1551 0.1501 0.1307 0.1202 0.0724 0.0637 0.0412 0.0315 0.0191 0.0114 0.0052 0.0120 0.0000 0.0000 0.0000 0.0000 0.0000 - 1980 1 1 1 0 4 0.0171 0.0527 0.1148 0.1286 0.1569 0.1335 0.1136 0.0812 0.0615 0.0435 0.0335 0.0215 0.0169 0.0097 0.0149 0.0000 0.0000 0.0000 0.0000 0.0000 - 1981 1 1 1 0 4 0.0220 0.0456 0.1138 0.1632 0.1384 0.1474 0.1121 0.0817 0.0566 0.0361 0.0272 0.0190 0.0139 0.0063 0.0168 0.0000 0.0000 0.0000 0.0000 0.0000 - 1982 1 1 1 0 4 0.0131 0.0480 0.0903 0.1633 0.1685 0.1323 0.1253 0.0874 0.0595 0.0370 0.0229 0.0197 0.0109 0.0080 0.0138 0.0000 0.0000 0.0000 0.0000 0.0000 - 1983 1 1 1 0 4 0.0094 0.0373 0.0972 0.1388 0.2068 0.1435 0.1209 0.0843 0.0544 0.0376 0.0263 0.0150 0.0102 0.0054 0.0128 0.0000 0.0000 0.0000 0.0000 0.0000 - 1984 1 1 1 0 4 0.0120 0.0614 0.1180 0.1511 0.1542 0.1752 0.1124 0.0727 0.0566 0.0305 0.0188 0.0143 0.0077 0.0046 0.0105 0.0000 0.0000 0.0000 0.0000 0.0000 - 1985 1 1 1 0 4 0.0064 0.0295 0.1191 0.1533 0.1683 0.1511 0.1451 0.0750 0.0484 0.0363 0.0257 0.0128 0.0099 0.0058 0.0131 0.0000 0.0000 0.0000 0.0000 0.0000 - 1986 1 1 1 0 4 0.0035 0.0274 0.0908 0.1844 0.1643 0.1523 0.1312 0.0982 0.0544 0.0326 0.0260 0.0119 0.0081 0.0062 0.0085 0.0000 0.0000 0.0000 0.0000 0.0000 - 1987 1 1 1 0 4 0.0038 0.0346 0.0878 0.1291 0.1987 0.1637 0.1203 0.0885 0.0678 0.0398 0.0239 0.0144 0.0100 0.0058 0.0117 0.0000 0.0000 0.0000 0.0000 0.0000 - 1988 1 1 1 0 4 0.0048 0.0317 0.1117 0.1593 0.1499 0.1979 0.1123 0.0701 0.0595 0.0443 0.0198 0.0132 0.0088 0.0054 0.0112 0.0000 0.0000 0.0000 0.0000 0.0000 - 1989 1 1 1 0 4 0.0077 0.0178 0.0686 0.1621 0.1713 0.1401 0.1707 0.0867 0.0584 0.0439 0.0315 0.0149 0.0108 0.0051 0.0102 0.0000 0.0000 0.0000 0.0000 0.0000 - 1990 1 1 1 0 4 0.0088 0.0221 0.0455 0.0972 0.1859 0.1656 0.1482 0.1348 0.0692 0.0449 0.0338 0.0184 0.0087 0.0060 0.0109 0.0000 0.0000 0.0000 0.0000 0.0000 - 1991 1 1 1 0 4 0.0045 0.0260 0.0574 0.0701 0.1216 0.1888 0.1699 0.1250 0.1056 0.0553 0.0292 0.0190 0.0126 0.0050 0.0101 0.0000 0.0000 0.0000 0.0000 0.0000 - 1992 1 1 1 0 4 0.0028 0.0162 0.0566 0.1095 0.1007 0.1461 0.1723 0.1291 0.0923 0.0753 0.0395 0.0227 0.0164 0.0077 0.0128 0.0000 0.0000 0.0000 0.0000 0.0000 - 1993 1 1 1 0 4 0.0035 0.0074 0.0348 0.0823 0.1371 0.1219 0.1291 0.1459 0.1197 0.0798 0.0592 0.0330 0.0172 0.0110 0.0183 0.0000 0.0000 0.0000 0.0000 0.0000 - 1994 1 1 1 0 4 0.0066 0.0133 0.0265 0.0490 0.0954 0.1461 0.1318 0.1284 0.1330 0.1006 0.0695 0.0439 0.0246 0.0109 0.0203 0.0000 0.0000 0.0000 0.0000 0.0000 - 1995 1 1 1 0 4 0.0059 0.0177 0.0599 0.0571 0.0885 0.1391 0.1623 0.1183 0.1125 0.0843 0.0538 0.0346 0.0259 0.0132 0.0268 0.0000 0.0000 0.0000 0.0000 0.0000 - 1996 1 1 1 0 4 0.0030 0.0126 0.0449 0.1225 0.0942 0.1136 0.1509 0.1414 0.0917 0.0690 0.0579 0.0397 0.0223 0.0134 0.0231 0.0000 0.0000 0.0000 0.0000 0.0000 - 1997 1 1 1 0 4 0.0052 0.0160 0.0423 0.1020 0.1781 0.1113 0.1204 0.1149 0.0965 0.0565 0.0464 0.0383 0.0236 0.0179 0.0306 0.0000 0.0000 0.0000 0.0000 0.0000 - 1998 1 1 1 0 4 0.0022 0.0088 0.0282 0.0544 0.1334 0.2181 0.1091 0.0981 0.0976 0.0722 0.0467 0.0349 0.0321 0.0206 0.0436 0.0000 0.0000 0.0000 0.0000 0.0000 - 1999 1 1 1 0 4 0.0026 0.0060 0.0170 0.0414 0.0691 0.1477 0.2165 0.0981 0.0842 0.0791 0.0655 0.0411 0.0380 0.0291 0.0647 0.0000 0.0000 0.0000 0.0000 0.0000 - 2000 1 1 1 0 4 0.0051 0.0075 0.0211 0.0338 0.0632 0.0930 0.1460 0.2367 0.0879 0.0580 0.0510 0.0501 0.0333 0.0263 0.0869 0.0000 0.0000 0.0000 0.0000 0.0000 - 2001 1 1 1 0 4 0.0018 0.0083 0.0173 0.0353 0.0488 0.0788 0.1057 0.1701 0.1959 0.0840 0.0458 0.0422 0.0363 0.0291 0.1006 0.0000 0.0000 0.0000 0.0000 0.0000 + 1935 1 1 1 0 4 0.1085 0.1645 0.1479 0.1896 0.1486 0.0762 0.0610 0.0375 0.0225 0.0165 0.0112 0.0055 0.0044 0.0021 0.0041 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1936 1 1 1 0 4 0.0729 0.1376 0.1219 0.1225 0.1496 0.1208 0.0864 0.0673 0.0388 0.0261 0.0183 0.0132 0.0099 0.0042 0.0105 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1937 1 1 1 0 4 0.0715 0.0929 0.1694 0.1244 0.1129 0.1208 0.1000 0.0802 0.0379 0.0364 0.0168 0.0108 0.0073 0.0060 0.0127 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1938 1 1 1 0 4 0.1129 0.1154 0.1649 0.1634 0.0741 0.0772 0.0835 0.0751 0.0559 0.0315 0.0184 0.0100 0.0049 0.0050 0.0080 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1939 1 1 1 0 4 0.0883 0.0863 0.1537 0.1421 0.1535 0.0809 0.0819 0.0690 0.0608 0.0323 0.0269 0.0116 0.0074 0.0010 0.0041 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1940 1 1 1 0 4 0.1354 0.0819 0.0997 0.1274 0.1618 0.1236 0.0779 0.0485 0.0484 0.0475 0.0202 0.0132 0.0081 0.0016 0.0049 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1941 1 1 1 0 4 0.1321 0.0840 0.1386 0.1137 0.1337 0.1125 0.1064 0.0596 0.0347 0.0352 0.0245 0.0127 0.0054 0.0022 0.0046 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1942 1 1 1 0 4 0.0924 0.1572 0.1021 0.1197 0.1161 0.1057 0.1015 0.0787 0.0431 0.0260 0.0257 0.0132 0.0100 0.0038 0.0047 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1943 1 1 1 0 4 0.0954 0.1839 0.2094 0.0969 0.0855 0.0744 0.0720 0.0754 0.0382 0.0258 0.0190 0.0104 0.0057 0.0036 0.0044 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1944 1 1 1 0 4 0.0616 0.0793 0.1749 0.1459 0.1004 0.0857 0.0946 0.0978 0.0637 0.0371 0.0230 0.0160 0.0071 0.0078 0.0051 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1945 1 1 1 0 4 0.0892 0.1572 0.1585 0.1452 0.0999 0.0767 0.0744 0.0780 0.0465 0.0273 0.0177 0.0127 0.0056 0.0055 0.0057 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1946 1 1 1 0 4 0.0457 0.1226 0.1913 0.1477 0.1363 0.0852 0.0767 0.0761 0.0457 0.0264 0.0177 0.0122 0.0056 0.0051 0.0057 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1947 1 1 1 0 4 0.0630 0.0716 0.1631 0.1503 0.1277 0.1112 0.0933 0.0852 0.0513 0.0297 0.0202 0.0150 0.0066 0.0058 0.0061 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1948 1 1 1 0 4 0.0630 0.0871 0.1647 0.1762 0.1370 0.0976 0.0798 0.0761 0.0457 0.0262 0.0181 0.0123 0.0055 0.0052 0.0054 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1949 1 1 1 0 4 0.0371 0.0577 0.1235 0.1477 0.1650 0.1137 0.1188 0.0934 0.0551 0.0264 0.0190 0.0173 0.0092 0.0058 0.0101 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1950 1 1 1 0 4 0.0190 0.0408 0.1129 0.1959 0.1516 0.1165 0.1108 0.0877 0.0771 0.0301 0.0239 0.0126 0.0063 0.0046 0.0101 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1951 1 1 1 0 4 0.0085 0.0519 0.1064 0.1568 0.1718 0.1291 0.1021 0.0955 0.0615 0.0518 0.0262 0.0169 0.0101 0.0041 0.0073 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1952 1 1 1 0 4 0.0226 0.0216 0.1210 0.1495 0.1397 0.1443 0.1017 0.0956 0.0746 0.0522 0.0331 0.0204 0.0094 0.0064 0.0079 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1953 1 1 1 0 4 0.0289 0.0584 0.0866 0.1703 0.1513 0.1355 0.1093 0.0797 0.0571 0.0451 0.0322 0.0181 0.0141 0.0039 0.0095 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1954 1 1 1 0 4 0.0387 0.0467 0.1034 0.1053 0.1854 0.1143 0.0980 0.0848 0.0646 0.0505 0.0429 0.0273 0.0173 0.0093 0.0116 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1955 1 1 1 0 4 0.0889 0.0744 0.0985 0.1275 0.1059 0.1523 0.1011 0.0768 0.0511 0.0404 0.0303 0.0231 0.0125 0.0079 0.0095 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1956 1 1 1 0 4 0.0940 0.1027 0.1119 0.0929 0.1271 0.0856 0.1181 0.0770 0.0622 0.0410 0.0257 0.0208 0.0151 0.0121 0.0139 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1957 1 1 1 0 4 0.0812 0.1069 0.1202 0.1188 0.0805 0.1076 0.0752 0.1030 0.0603 0.0487 0.0323 0.0228 0.0172 0.0097 0.0155 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1958 1 1 1 0 4 0.1450 0.1145 0.1120 0.1036 0.0985 0.0761 0.0772 0.0571 0.0762 0.0471 0.0288 0.0209 0.0163 0.0111 0.0154 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1959 1 1 1 0 4 0.1324 0.0669 0.1993 0.1381 0.1077 0.0862 0.0522 0.0593 0.0388 0.0387 0.0234 0.0180 0.0115 0.0105 0.0168 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1960 1 1 1 0 4 0.1126 0.0778 0.1170 0.2314 0.1287 0.0831 0.0570 0.0437 0.0374 0.0305 0.0256 0.0173 0.0120 0.0086 0.0171 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1961 1 1 1 0 4 0.0644 0.1012 0.1132 0.1301 0.2120 0.1110 0.0719 0.0492 0.0351 0.0305 0.0239 0.0200 0.0124 0.0087 0.0164 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1962 1 1 1 0 4 0.0655 0.0763 0.1229 0.1351 0.1461 0.1892 0.0738 0.0491 0.0333 0.0261 0.0252 0.0192 0.0131 0.0080 0.0172 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1963 1 1 1 0 4 0.0818 0.0735 0.1434 0.1511 0.1308 0.1204 0.1309 0.0533 0.0279 0.0179 0.0172 0.0150 0.0121 0.0086 0.0160 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1964 1 1 1 0 4 0.0504 0.0929 0.1002 0.1885 0.1384 0.1168 0.0987 0.0976 0.0381 0.0189 0.0146 0.0121 0.0086 0.0080 0.0161 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1965 1 1 1 0 4 0.0524 0.0817 0.1552 0.1160 0.1898 0.1161 0.0880 0.0725 0.0562 0.0209 0.0129 0.0093 0.0083 0.0058 0.0148 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1966 1 1 1 0 4 0.0471 0.0494 0.1096 0.1534 0.1143 0.1619 0.1010 0.0805 0.0661 0.0464 0.0217 0.0137 0.0106 0.0056 0.0187 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1967 1 1 1 0 4 0.0778 0.0711 0.0881 0.1643 0.1548 0.0963 0.1196 0.0687 0.0479 0.0374 0.0269 0.0136 0.0094 0.0067 0.0173 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1968 1 1 1 0 4 0.0777 0.1108 0.0990 0.1021 0.1611 0.1256 0.0794 0.0932 0.0489 0.0354 0.0268 0.0173 0.0080 0.0044 0.0102 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1969 1 1 1 0 4 0.1105 0.0591 0.1727 0.0987 0.0872 0.1288 0.1020 0.0659 0.0637 0.0329 0.0272 0.0203 0.0117 0.0066 0.0128 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1970 1 1 1 0 4 0.1088 0.0799 0.0987 0.1886 0.0915 0.0820 0.1064 0.0727 0.0506 0.0450 0.0242 0.0158 0.0122 0.0089 0.0147 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1971 1 1 1 0 4 0.1255 0.1017 0.1335 0.1164 0.1729 0.0706 0.0621 0.0688 0.0479 0.0387 0.0238 0.0147 0.0094 0.0061 0.0081 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1972 1 1 1 0 4 0.1139 0.1099 0.1348 0.1347 0.1125 0.1493 0.0582 0.0573 0.0440 0.0318 0.0209 0.0132 0.0065 0.0049 0.0082 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1973 1 1 1 0 4 0.0269 0.0661 0.1138 0.1376 0.1489 0.1187 0.1492 0.0621 0.0486 0.0480 0.0308 0.0160 0.0136 0.0054 0.0144 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1974 1 1 1 0 4 0.0202 0.0454 0.1028 0.1246 0.1325 0.1304 0.1035 0.1276 0.0472 0.0465 0.0485 0.0268 0.0151 0.0115 0.0173 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1975 1 1 1 0 4 0.0178 0.0423 0.1025 0.1322 0.1423 0.1332 0.1220 0.0778 0.0906 0.0343 0.0321 0.0265 0.0170 0.0080 0.0215 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1976 1 1 1 0 4 0.0491 0.0539 0.1244 0.1442 0.1340 0.1095 0.1065 0.0744 0.0602 0.0617 0.0210 0.0196 0.0148 0.0081 0.0185 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1977 1 1 1 0 4 0.0288 0.0583 0.1065 0.1362 0.1577 0.1079 0.1049 0.0848 0.0654 0.0491 0.0453 0.0172 0.0137 0.0075 0.0165 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1978 1 1 1 0 4 0.0390 0.0691 0.1275 0.1387 0.1465 0.1236 0.1018 0.0691 0.0590 0.0425 0.0320 0.0190 0.0103 0.0065 0.0153 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1979 1 1 1 0 4 0.0165 0.0644 0.1065 0.1551 0.1501 0.1307 0.1202 0.0724 0.0637 0.0412 0.0315 0.0191 0.0114 0.0052 0.0120 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1980 1 1 1 0 4 0.0171 0.0527 0.1148 0.1286 0.1569 0.1335 0.1136 0.0812 0.0615 0.0435 0.0335 0.0215 0.0169 0.0097 0.0149 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1981 1 1 1 0 4 0.0220 0.0456 0.1138 0.1632 0.1384 0.1474 0.1121 0.0817 0.0566 0.0361 0.0272 0.0190 0.0139 0.0063 0.0168 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1982 1 1 1 0 4 0.0131 0.0480 0.0903 0.1633 0.1685 0.1323 0.1253 0.0874 0.0595 0.0370 0.0229 0.0197 0.0109 0.0080 0.0138 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1983 1 1 1 0 4 0.0094 0.0373 0.0972 0.1388 0.2068 0.1435 0.1209 0.0843 0.0544 0.0376 0.0263 0.0150 0.0102 0.0054 0.0128 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1984 1 1 1 0 4 0.0120 0.0614 0.1180 0.1511 0.1542 0.1752 0.1124 0.0727 0.0566 0.0305 0.0188 0.0143 0.0077 0.0046 0.0105 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1985 1 1 1 0 4 0.0064 0.0295 0.1191 0.1533 0.1683 0.1511 0.1451 0.0750 0.0484 0.0363 0.0257 0.0128 0.0099 0.0058 0.0131 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1986 1 1 1 0 4 0.0035 0.0274 0.0908 0.1844 0.1643 0.1523 0.1312 0.0982 0.0544 0.0326 0.0260 0.0119 0.0081 0.0062 0.0085 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1987 1 1 1 0 4 0.0038 0.0346 0.0878 0.1291 0.1987 0.1637 0.1203 0.0885 0.0678 0.0398 0.0239 0.0144 0.0100 0.0058 0.0117 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1988 1 1 1 0 4 0.0048 0.0317 0.1117 0.1593 0.1499 0.1979 0.1123 0.0701 0.0595 0.0443 0.0198 0.0132 0.0088 0.0054 0.0112 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1989 1 1 1 0 4 0.0077 0.0178 0.0686 0.1621 0.1713 0.1401 0.1707 0.0867 0.0584 0.0439 0.0315 0.0149 0.0108 0.0051 0.0102 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1990 1 1 1 0 4 0.0088 0.0221 0.0455 0.0972 0.1859 0.1656 0.1482 0.1348 0.0692 0.0449 0.0338 0.0184 0.0087 0.0060 0.0109 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1991 1 1 1 0 4 0.0045 0.0260 0.0574 0.0701 0.1216 0.1888 0.1699 0.1250 0.1056 0.0553 0.0292 0.0190 0.0126 0.0050 0.0101 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1992 1 1 1 0 4 0.0028 0.0162 0.0566 0.1095 0.1007 0.1461 0.1723 0.1291 0.0923 0.0753 0.0395 0.0227 0.0164 0.0077 0.0128 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1993 1 1 1 0 4 0.0035 0.0074 0.0348 0.0823 0.1371 0.1219 0.1291 0.1459 0.1197 0.0798 0.0592 0.0330 0.0172 0.0110 0.0183 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1994 1 1 1 0 4 0.0066 0.0133 0.0265 0.0490 0.0954 0.1461 0.1318 0.1284 0.1330 0.1006 0.0695 0.0439 0.0246 0.0109 0.0203 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1995 1 1 1 0 4 0.0059 0.0177 0.0599 0.0571 0.0885 0.1391 0.1623 0.1183 0.1125 0.0843 0.0538 0.0346 0.0259 0.0132 0.0268 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1996 1 1 1 0 4 0.0030 0.0126 0.0449 0.1225 0.0942 0.1136 0.1509 0.1414 0.0917 0.0690 0.0579 0.0397 0.0223 0.0134 0.0231 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1997 1 1 1 0 4 0.0052 0.0160 0.0423 0.1020 0.1781 0.1113 0.1204 0.1149 0.0965 0.0565 0.0464 0.0383 0.0236 0.0179 0.0306 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1998 1 1 1 0 4 0.0022 0.0088 0.0282 0.0544 0.1334 0.2181 0.1091 0.0981 0.0976 0.0722 0.0467 0.0349 0.0321 0.0206 0.0436 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1999 1 1 1 0 4 0.0026 0.0060 0.0170 0.0414 0.0691 0.1477 0.2165 0.0981 0.0842 0.0791 0.0655 0.0411 0.0380 0.0291 0.0647 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2000 1 1 1 0 4 0.0051 0.0075 0.0211 0.0338 0.0632 0.0930 0.1460 0.2367 0.0879 0.0580 0.0510 0.0501 0.0333 0.0263 0.0869 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2001 1 1 1 0 4 0.0018 0.0083 0.0173 0.0353 0.0488 0.0788 0.1057 0.1701 0.1959 0.0840 0.0458 0.0422 0.0363 0.0291 0.1006 ## 0.0000 0.0000 0.0000 0.0000 0.0000 2002 1 1 1 0 2 0.0038 0.0129 0.0380 0.0326 0.0506 0.0570 0.0721 0.0877 0.1285 0.1691 0.0905 0.0510 0.0453 0.0401 0.0278 0.0217 0.0190 0.0163 0.0111 0.0250 2003 1 1 1 0 2 0.0058 0.0160 0.0477 0.0811 0.0562 0.0595 0.0680 0.0688 0.0887 0.1136 0.1284 0.0758 0.0423 0.0350 0.0341 0.0208 0.0119 0.0141 0.0102 0.0220 2004 1 1 1 0 2 0.0158 0.0156 0.0371 0.0902 0.1056 0.0659 0.0641 0.0679 0.0615 0.0782 0.0985 0.1008 0.0536 0.0346 0.0276 0.0280 0.0161 0.0099 0.0091 0.0200 @@ -616,19 +616,19 @@ 2014 1 1 1 0 2 0.0022 0.0043 0.0144 0.0582 0.1035 0.1241 0.1569 0.1317 0.1078 0.0915 0.0529 0.0354 0.0285 0.0174 0.0149 0.0086 0.0081 0.0073 0.0071 0.0252 ## survey female ## year gear area group sex age_err | data columns (numbers or proportions) - 1980 6 1 1 1 4 0.0292 0.0363 0.0740 0.0861 0.0917 0.0755 0.0688 0.0551 0.0454 0.0225 0.0187 0.0188 0.0090 0.0065 0.0150 0.0000 0.0000 0.0000 0.0000 0.0000 - 1981 6 1 1 1 4 0.0535 0.0246 0.0522 0.0731 0.0711 0.0634 0.0664 0.0621 0.0432 0.0307 0.0246 0.0198 0.0144 0.0082 0.0151 0.0000 0.0000 0.0000 0.0000 0.0000 - 1982 6 1 1 1 4 0.0462 0.0467 0.0630 0.1072 0.0994 0.0792 0.0706 0.0509 0.0279 0.0177 0.0121 0.0078 0.0076 0.0080 0.0121 0.0000 0.0000 0.0000 0.0000 0.0000 - 1983 6 1 1 1 4 0.0584 0.0556 0.0608 0.0604 0.0723 0.0690 0.0580 0.0665 0.0407 0.0298 0.0157 0.0102 0.0078 0.0051 0.0119 0.0000 0.0000 0.0000 0.0000 0.0000 - 1984 6 1 1 1 4 0.0295 0.0444 0.0566 0.0505 0.0521 0.0720 0.0505 0.0319 0.0320 0.0229 0.0143 0.0123 0.0091 0.0088 0.0196 0.0000 0.0000 0.0000 0.0000 0.0000 - 1985 6 1 1 1 4 0.0308 0.0280 0.0734 0.0625 0.0594 0.0459 0.0549 0.0457 0.0316 0.0257 0.0108 0.0095 0.0049 0.0062 0.0218 0.0000 0.0000 0.0000 0.0000 0.0000 - 1986 6 1 1 1 4 0.0267 0.0306 0.0375 0.0589 0.0553 0.0510 0.0461 0.0539 0.0402 0.0333 0.0227 0.0140 0.0121 0.0086 0.0380 0.0000 0.0000 0.0000 0.0000 0.0000 - 1996 6 1 1 1 4 0.0126 0.0207 0.0567 0.0961 0.0524 0.0694 0.0836 0.0665 0.0508 0.0342 0.0261 0.0191 0.0090 0.0065 0.0077 0.0000 0.0000 0.0000 0.0000 0.0000 - 1997 6 1 1 1 4 0.0057 0.0173 0.0262 0.0646 0.1270 0.0654 0.0628 0.0581 0.0481 0.0232 0.0212 0.0169 0.0096 0.0082 0.0076 0.0000 0.0000 0.0000 0.0000 0.0000 + 1980 6 1 1 1 4 0.0292 0.0363 0.0740 0.0861 0.0917 0.0755 0.0688 0.0551 0.0454 0.0225 0.0187 0.0188 0.0090 0.0065 0.0150 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1981 6 1 1 1 4 0.0535 0.0246 0.0522 0.0731 0.0711 0.0634 0.0664 0.0621 0.0432 0.0307 0.0246 0.0198 0.0144 0.0082 0.0151 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1982 6 1 1 1 4 0.0462 0.0467 0.0630 0.1072 0.0994 0.0792 0.0706 0.0509 0.0279 0.0177 0.0121 0.0078 0.0076 0.0080 0.0121 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1983 6 1 1 1 4 0.0584 0.0556 0.0608 0.0604 0.0723 0.0690 0.0580 0.0665 0.0407 0.0298 0.0157 0.0102 0.0078 0.0051 0.0119 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1984 6 1 1 1 4 0.0295 0.0444 0.0566 0.0505 0.0521 0.0720 0.0505 0.0319 0.0320 0.0229 0.0143 0.0123 0.0091 0.0088 0.0196 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1985 6 1 1 1 4 0.0308 0.0280 0.0734 0.0625 0.0594 0.0459 0.0549 0.0457 0.0316 0.0257 0.0108 0.0095 0.0049 0.0062 0.0218 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1986 6 1 1 1 4 0.0267 0.0306 0.0375 0.0589 0.0553 0.0510 0.0461 0.0539 0.0402 0.0333 0.0227 0.0140 0.0121 0.0086 0.0380 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1996 6 1 1 1 4 0.0126 0.0207 0.0567 0.0961 0.0524 0.0694 0.0836 0.0665 0.0508 0.0342 0.0261 0.0191 0.0090 0.0065 0.0077 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1997 6 1 1 1 4 0.0057 0.0173 0.0262 0.0646 0.1270 0.0654 0.0628 0.0581 0.0481 0.0232 0.0212 0.0169 0.0096 0.0082 0.0076 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1999 6 1 1 1 4 0.0068 0.0144 0.0228 0.0390 0.0491 0.0886 0.1125 0.0564 0.0476 0.0332 0.0221 0.0130 0.0102 0.0067 0.0117 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2000 6 1 1 1 4 0.0118 0.0110 0.0226 0.0361 0.0465 0.0594 0.0897 0.1106 0.0455 0.0308 0.0248 0.0189 0.0115 0.0068 0.0175 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2001 6 1 1 1 4 0.0096 0.0166 0.0180 0.0248 0.0339 0.0519 0.0549 0.0833 0.0921 0.0572 0.0333 0.0187 0.0136 0.0087 0.0118 ## 0.0000 0.0000 0.0000 0.0000 0.0000 1998 6 1 1 1 2 0.0069 0.0121 0.0177 0.0272 0.0487 0.1145 0.0622 0.0447 0.0528 0.0368 0.0225 0.0153 0.0161 0.0137 0.0095 0.0053 0.0037 0.0015 0.0026 0.0043 - 1999 6 1 1 1 4 0.0068 0.0144 0.0228 0.0390 0.0491 0.0886 0.1125 0.0564 0.0476 0.0332 0.0221 0.0130 0.0102 0.0067 0.0117 0.0000 0.0000 0.0000 0.0000 0.0000 - 2000 6 1 1 1 4 0.0118 0.0110 0.0226 0.0361 0.0465 0.0594 0.0897 0.1106 0.0455 0.0308 0.0248 0.0189 0.0115 0.0068 0.0175 0.0000 0.0000 0.0000 0.0000 0.0000 - 2001 6 1 1 1 4 0.0096 0.0166 0.0180 0.0248 0.0339 0.0519 0.0549 0.0833 0.0921 0.0572 0.0333 0.0187 0.0136 0.0087 0.0118 0.0000 0.0000 0.0000 0.0000 0.0000 2002 6 1 1 1 2 0.0113 0.0303 0.0418 0.0265 0.0319 0.0358 0.0443 0.0440 0.0666 0.0760 0.0447 0.0259 0.0155 0.0139 0.0103 0.0062 0.0063 0.0042 0.0030 0.0058 2003 6 1 1 1 2 0.0170 0.0286 0.0523 0.0662 0.0394 0.0300 0.0392 0.0364 0.0440 0.0501 0.0534 0.0314 0.0170 0.0119 0.0088 0.0085 0.0036 0.0045 0.0031 0.0057 2004 6 1 1 1 2 0.0233 0.0264 0.0438 0.0678 0.0555 0.0324 0.0339 0.0275 0.0295 0.0339 0.0334 0.0299 0.0161 0.0097 0.0064 0.0055 0.0033 0.0023 0.0016 0.0039 @@ -644,19 +644,19 @@ 2014 6 1 1 1 2 0.0051 0.0119 0.0409 0.0986 0.0995 0.0750 0.0742 0.0565 0.0488 0.0339 0.0141 0.0098 0.0065 0.0047 0.0026 0.0014 0.0008 0.0010 0.0005 0.0030 ## survey male ## year gear area group sex age_err | data columns (numbers or proportions) - 1980 6 1 1 2 4 0.0174 0.0433 0.0462 0.0488 0.0444 0.0482 0.0315 0.0200 0.0155 0.0091 0.0077 0.0031 0.0036 0.0037 0.0048 0.0000 0.0000 0.0000 0.0000 0.0000 - 1981 6 1 1 2 4 0.0503 0.0396 0.0576 0.0501 0.0457 0.0435 0.0305 0.0225 0.0150 0.0072 0.0075 0.0013 0.0028 0.0018 0.0022 0.0000 0.0000 0.0000 0.0000 0.0000 - 1982 6 1 1 2 4 0.0270 0.0313 0.0439 0.0722 0.0461 0.0314 0.0308 0.0210 0.0149 0.0105 0.0064 0.0028 0.0020 0.0011 0.0022 0.0000 0.0000 0.0000 0.0000 0.0000 - 1983 6 1 1 2 4 0.0369 0.0408 0.0429 0.0502 0.0634 0.0496 0.0274 0.0241 0.0123 0.0096 0.0091 0.0041 0.0021 0.0023 0.0031 0.0000 0.0000 0.0000 0.0000 0.0000 - 1984 6 1 1 2 4 0.0399 0.0566 0.0553 0.0642 0.0662 0.0692 0.0458 0.0286 0.0240 0.0163 0.0106 0.0064 0.0020 0.0022 0.0064 0.0000 0.0000 0.0000 0.0000 0.0000 - 1985 6 1 1 2 4 0.0241 0.0415 0.0811 0.0669 0.0609 0.0574 0.0644 0.0387 0.0188 0.0113 0.0089 0.0067 0.0030 0.0033 0.0023 0.0000 0.0000 0.0000 0.0000 0.0000 - 1986 6 1 1 2 4 0.0283 0.0281 0.0394 0.0761 0.0503 0.0452 0.0542 0.0479 0.0334 0.0216 0.0153 0.0121 0.0070 0.0041 0.0082 0.0000 0.0000 0.0000 0.0000 0.0000 - 1996 6 1 1 2 4 0.0062 0.0096 0.0277 0.0510 0.0323 0.0390 0.0490 0.0426 0.0339 0.0303 0.0244 0.0178 0.0110 0.0075 0.0061 0.0000 0.0000 0.0000 0.0000 0.0000 - 1997 6 1 1 2 4 0.0035 0.0073 0.0178 0.0400 0.0871 0.0390 0.0436 0.0446 0.0437 0.0290 0.0267 0.0207 0.0133 0.0076 0.0141 0.0000 0.0000 0.0000 0.0000 0.0000 + 1980 6 1 1 2 4 0.0174 0.0433 0.0462 0.0488 0.0444 0.0482 0.0315 0.0200 0.0155 0.0091 0.0077 0.0031 0.0036 0.0037 0.0048 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1981 6 1 1 2 4 0.0503 0.0396 0.0576 0.0501 0.0457 0.0435 0.0305 0.0225 0.0150 0.0072 0.0075 0.0013 0.0028 0.0018 0.0022 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1982 6 1 1 2 4 0.0270 0.0313 0.0439 0.0722 0.0461 0.0314 0.0308 0.0210 0.0149 0.0105 0.0064 0.0028 0.0020 0.0011 0.0022 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1983 6 1 1 2 4 0.0369 0.0408 0.0429 0.0502 0.0634 0.0496 0.0274 0.0241 0.0123 0.0096 0.0091 0.0041 0.0021 0.0023 0.0031 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1984 6 1 1 2 4 0.0399 0.0566 0.0553 0.0642 0.0662 0.0692 0.0458 0.0286 0.0240 0.0163 0.0106 0.0064 0.0020 0.0022 0.0064 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1985 6 1 1 2 4 0.0241 0.0415 0.0811 0.0669 0.0609 0.0574 0.0644 0.0387 0.0188 0.0113 0.0089 0.0067 0.0030 0.0033 0.0023 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1986 6 1 1 2 4 0.0283 0.0281 0.0394 0.0761 0.0503 0.0452 0.0542 0.0479 0.0334 0.0216 0.0153 0.0121 0.0070 0.0041 0.0082 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1996 6 1 1 2 4 0.0062 0.0096 0.0277 0.0510 0.0323 0.0390 0.0490 0.0426 0.0339 0.0303 0.0244 0.0178 0.0110 0.0075 0.0061 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1997 6 1 1 2 4 0.0035 0.0073 0.0178 0.0400 0.0871 0.0390 0.0436 0.0446 0.0437 0.0290 0.0267 0.0207 0.0133 0.0076 0.0141 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 1999 6 1 1 2 4 0.0041 0.0070 0.0126 0.0247 0.0383 0.0742 0.1031 0.0452 0.0380 0.0332 0.0295 0.0180 0.0134 0.0115 0.0132 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2000 6 1 1 2 4 0.0055 0.0070 0.0155 0.0197 0.0345 0.0475 0.0747 0.0836 0.0375 0.0316 0.0323 0.0219 0.0169 0.0102 0.0184 ## 0.0000 0.0000 0.0000 0.0000 0.0000 + 2001 6 1 1 2 4 0.0067 0.0093 0.0104 0.0171 0.0243 0.0330 0.0548 0.0764 0.0774 0.0460 0.0386 0.0294 0.0209 0.0135 0.0137 ## 0.0000 0.0000 0.0000 0.0000 0.0000 1998 6 1 1 2 2 0.0030 0.0035 0.0076 0.0170 0.0353 0.0716 0.0437 0.0319 0.0388 0.0377 0.0302 0.0323 0.0320 0.0261 0.0242 0.0167 0.0077 0.0058 0.0057 0.0111 - 1999 6 1 1 2 4 0.0041 0.0070 0.0126 0.0247 0.0383 0.0742 0.1031 0.0452 0.0380 0.0332 0.0295 0.0180 0.0134 0.0115 0.0132 0.0000 0.0000 0.0000 0.0000 0.0000 - 2000 6 1 1 2 4 0.0055 0.0070 0.0155 0.0197 0.0345 0.0475 0.0747 0.0836 0.0375 0.0316 0.0323 0.0219 0.0169 0.0102 0.0184 0.0000 0.0000 0.0000 0.0000 0.0000 - 2001 6 1 1 2 4 0.0067 0.0093 0.0104 0.0171 0.0243 0.0330 0.0548 0.0764 0.0774 0.0460 0.0386 0.0294 0.0209 0.0135 0.0137 0.0000 0.0000 0.0000 0.0000 0.0000 2002 6 1 1 2 2 0.0067 0.0139 0.0249 0.0166 0.0220 0.0222 0.0285 0.0326 0.0463 0.0709 0.0368 0.0244 0.0222 0.0233 0.0143 0.0112 0.0107 0.0095 0.0058 0.0129 2003 6 1 1 2 2 0.0062 0.0147 0.0324 0.0308 0.0227 0.0205 0.0266 0.0252 0.0332 0.0456 0.0489 0.0342 0.0209 0.0160 0.0183 0.0130 0.0109 0.0066 0.0066 0.0159 2004 6 1 1 2 2 0.0122 0.0134 0.0279 0.0516 0.0442 0.0224 0.0251 0.0303 0.0339 0.0358 0.0483 0.0510 0.0313 0.0185 0.0172 0.0138 0.0115 0.0057 0.0059 0.0138 diff --git a/fba/HALIBUT/R/R4iSCAM.R b/fba/HALIBUT/R/R4iSCAM.R index 6d8fd4ad..7205fdc9 100644 --- a/fba/HALIBUT/R/R4iSCAM.R +++ b/fba/HALIBUT/R/R4iSCAM.R @@ -27,6 +27,12 @@ setwd(.PWD) .THEME <- theme_bw(12) .UNITS <- "(mlb)" +# | Labels for gear sex area and group index. +.GEAR = c("Directed","Wasteage","Bycatch","Sport","Personal","Setline Survey") +.SEX = c("F & M","Female","Male") +.AREA = c("Coast wide") +.GROUP = c("Pacific Halibut") + .srcRlib <- function() { @@ -59,9 +65,6 @@ for(nm in .RFILES) source(file.path(.LIB, nm), echo=FALSE) .plotCatch( M ) .plotIndex( M ) .plotWeightAtAge( M ) -.plotAgeComps( M ) -.plotAgeCompResiduals( M ) -.plotAgeSummary( M ) .plotSpawnBiomass( M ) .plotDepletion( M ) .plotMortality( M ) @@ -69,9 +72,15 @@ for(nm in .RFILES) source(file.path(.LIB, nm), echo=FALSE) .plotStockRecruit( M ) .plotRecruitsPerSpawner( M ) .plotSurveyFit( M ) +.plotQ( M ) .plotCatchResidual( M ) .plotIndexResidual( M ) .plotRecruitmentResidual( M ) +.plotAgeComps( M ) +.plotAgeCompResiduals( M ) +.plotAgeSummary( M ) +.plotAgeBars( M ) +.plotSlx( M ) # |----------------------------------------------------------------------------------| # | guiView: Main function for starting the iSCAM Gui diff --git a/scripts/Template.ctl b/scripts/Template.ctl index bb28e38d..04b1801d 100644 --- a/scripts/Template.ctl +++ b/scripts/Template.ctl @@ -1,32 +1,40 @@ -## ------------------------------------------------------------------------- ## -## CONTROL FILE TEMPLATE ## -## ------------------------------------------------------------------------- ## -## -## -## ------------------------------------------------------------------------- ## -## CONTROLS FOR LEADING PARAMETERS ## -## Prior descriptions: ## -## -0 uniform (0,0) ## -## -1 normal (p1=mu,p2=sig) ## -## -2 lognormal (p1=log(mu),p2=sig) ## -## -3 beta (p1=alpha,p2=beta) ## -## -4 gamma (p1=alpha,p2=beta) ## -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## CONTROL FILE TEMPLATE ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + +## ———————————————————————————————————————————————————————————————————————————————————— ## +## CONTROLS FOR LEADING PARAMETERS ## +## Prior descriptions: ## +## -0 uniform (0,0) ## +## -1 normal (p1=mu,p2=sig) ## +## -2 lognormal (p1=log(mu),p2=sig) ## +## -3 beta (p1=alpha,p2=beta) ## +## -4 gamma (p1=alpha,p2=beta) ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## npar -7 -## ival lb ub phz prior p1 p2 #parameter ## - 1.0 0.0 10 2 0 0.0 10.0 #log_ro ## - 0.75 0.2 1.0 2 3 3.00 2.00 #steepness ## - -1.203973 -3.0 2.0 4 1 -1.203 0.15 #log_m g&b ## - 1.0 0.0 10 1 0 0.0 10 #log_avgrec ## - 1.0 0.0 10 1 0 0.0 10 #log_recinit ## - 0.5 0.01 0.99 3 3 3.00 5.00 #rho ## - 0.8 0.01 5.0 3 4 1.01 1.01 #vartheta ## -## ------------------------------------------------------------------------- ## -## -## ------------------------------------------------------------------------- ## -## CONTROL PARAMETERS FOR AGE/SIZE COMPOSITION DATA FOR na_gears ## -## ------------------------------------------------------------------------- ## + 7 +## ival lb ub phz prior p1 p2 #parameter ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + 1.0 0.0 10 2 0 0.0 10.0 #log_ro ## + 0.75 0.2 1.0 2 3 3.00 2.00 #steepness ## + -1.203973 -3.0 2.0 4 1 -1.203 0.15 #log_m g&b ## + 1.0 0.0 10 1 0 0.0 10 #log_avgrec ## + 1.0 0.0 10 1 0 0.0 10 #log_recinit ## + 0.00 0.00 1.00 -3 3 3.00 5.00 #rho ## + 0.6 0.01 5.0 -3 4 1.01 1.01 #sigma_r ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + +## ———————————————————————————————————————————————————————————————————————————————————— ## +## CONTROL PARAMETERS FOR AGE/SIZE COMPOSITION DATA FOR na_gears ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## Likelihood type for each gear: ## -1 : multivariate logistic (dmvlogistic) ## -2 : multinomial, sample size based on input data @@ -34,20 +42,76 @@ ## -4 : logistic_normal, AR1 ## -5 : logistic_normal, AR2 ## -6 : multinomial with estimated effective samples size (log_age_tau2 phz >0) -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## ## Number of columns == na_gears. - 1 ## Gear Index - 3 ## Likelihood type - 0.000 ## Minimum proportion for aggregation & tail compression - 0.0000 ## Small constant to add to comps & renormalize - 1 ## phase for log_age_tau2 estimation. - 2 ## phase for phi1 estimation: bounded (-1,1) AR1 - -2 ## phase for phi2 estimation: bounded (0,1) AR2 - -2 ## phase for degrees of freedom for student T. - -12345 ## int check (-12345) -## ------------------------------------------------------------------------- ## - + 1 ## • Gear Index + 3 ## • Likelihood type + 0.000 ## • Minimum proportion for aggregation & tail compression + 0.0000 ## • Small constant to add to comps & renormalize + 1 ## • phase for log_age_tau2 estimation. + 2 ## • phase for phi1 estimation: bounded (-1,1) AR1 + -2 ## • phase for phi2 estimation: bounded (0,1) AR2 + -2 ## • phase for degrees of freedom for student T. + -12345 ## • int check (-12345) +## ———————————————————————————————————————————————————————————————————————————————————— ## ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## SELECTIVITY CONTROLS ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## - Each gear must have at least one selectivity and retention curve. +## - Use a -ve phase with the gear index to mirror another gear. Note +## that if you mirror another gear, it must have the same sel type and +## age and year nodes so that the arrays are the same shape & block years. +## +## • Index = gear index for selectivity curve. +## • sel_type = type of selectivity function (see Legend). +## • sel_mu = mean age/length 50% selectivity. +## • sel_sd = std in 50% selectivity +## • sex_dep = 0 -> no; 1 -> offset for sex 2. +## • size_nodes = # of nodes for age/size cubic spline. +## • year_nodes = # of nodes for time varying bicubic spline. +## • phz_mirror = phase of estimation (-ve phase to mirror selextivity index) +## • lam1 = penalty weight for 2nd differences (w = 1/(2•sig^2)). +## • lam2 = penalty weight for dome-shaped selectivity. +## • lam3 = penalty weight for time-varying selectivity. +## • start_block = year index for first year of selectivity curve. +## ———————————————————————————————————————————————————————————————————————————————————— ## +## sel_nBlocks ret_nBlocks ## Gear Description. + 1 1 ## Commercial retained +## ———————————————————————————————————————————————————————————————————————————————————— ## +## Selectivity P(capture of all size/age) +## slx_dControls +## • index for sex (0=both, 1=female, 2=male) +## sel sel sel age year phz or start end ## +## Index type mu sd sex nodes nodes mirror lam1 lam2 lam3 | block block ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + 1 2 3.5 2.0 0 0 0 2 2.0 2.0 12.5 1968 1979 +## ———————————————————————————————————————————————————————————————————————————————————— ## +## Retention P(retaining size/age) +## ret_dControls +## • index for sex (0=both, 1=female, 2=male) +## sel sel sel age year phz or start end ## +## Index type mu sd sex nodes nodes mirror lam1 lam2 lam3 | block block ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + -1 2 10 2.0 1 0 0 -2 0.0 0.0 0.0 1968 1979 +## ———————————————————————————————————————————————————————————————————————————————————— ## +## LEGEND FOR SELECTIVITY TYPES (sel_type) ## +## ———————————————————————————————————————————————————————————————————————————————————— ## +## sel | No. | +## type | parameters | Description +## ———————————————————————————————————————————————————————————————————————————————————— ## +## 1 | 2 | • Logistic curve with mean and standard dev at age = p(50%) +## 2 | (nages-1) | • Age-specific selectivity coefficients for (sage):(nage-1) +## 3 | age_nodes | • Age-specific coefficients based on cubic-spline interpolation +## 4 | n*age_nodes| • Annual age-specific coeffs using cubic-spline interpolation +## 5 | nyr*nage | • Bicubic spline interpolation over time & age. +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + +## TO BE DEPRECATED ## ------------------------------------------------------------------------- ## ## SELECTIVITY PARAMETERS Columns for gear ## ## OPTIONS FOR SELECTIVITY (isel_type): ## @@ -64,16 +128,16 @@ ## sig=0.05 0.10 0.15 0.20 0.30 0.40 0.50 ## ## wt =200. 50.0 22.2 12.5 5.56 3.12 2.00 ## ## ------------------------------------------------------------------------- ## - 2 # 1 -selectivity type ivector(isel_type) for gear - 3.5 # 2 -Age/length at 50% selectivity (logistic) - 0.5 # 3 -STD at 50% selectivity (logistic) - 7 # 4 -No. of age nodes for each gear (0=ignore) - 12 # 5 -No. of year nodes for 2d spline(0=ignore) - 3 # 6 -Phase of estimation (-1 for fixed) - 2.00 # 7 -Penalty wt for 2nd differences w=1/(2*sig^2) - 2.00 # 8 -Penalty wt for dome-shaped w=1/(2*sig^2) - 12.5 # 9 -Penalty wt for time-varying selectivity - 1 # 10 -n_sel_blocks (number of selex blocks) + 2 # 1 -selectivity type ivector(isel_type) for gear + 3.5 # 2 -Age/length at 50% selectivity (logistic) + 0.5 # 3 -STD at 50% selectivity (logistic) + 7 # 4 -No. of age nodes for each gear (0=ignore) + 12 # 5 -No. of year nodes for 2d spline(0=ignore) + 3 # 6 -Phase of estimation (-1 for fixed) + 2.00 # 7 -Penalty wt for 2nd differences w=1/(2*sig^2) + 2.00 # 8 -Penalty wt for dome-shaped w=1/(2*sig^2) + 12.5 # 9 -Penalty wt for time-varying selectivity + 1 # 10 -n_sel_blocks (number of selex blocks) ## ------------------------------------------------------------------------- ## ## Start year of each time block: 1 row for each gear 1968 @@ -98,25 +162,33 @@ ## Year position of the knots (vector must be equal to the number of nodes) ## ———————————————————————————————————————————————————————————————————————————————————— ## -## -## -## ------------------------------------------------------------------------- ## -## PRIORS FOR SURVEY Q ## -## Prior type: ## -## 0 - uninformative prior ## -## 1 - normal prior density for log(q) ## -## 2 - random walk in q ## -## ------------------------------------------------------------------------- ## -1 # -number of surveys (nits) -0 # -prior type (see legend above) -0 # -prior log(mean) -0 # -prior sd -## ------------------------------------------------------------------------- ## -## -## ------------------------------------------------------------------------- ## -## OTHER MISCELANEOUS CONTROLS ## -## ------------------------------------------------------------------------- ## + + + + +## ———————————————————————————————————————————————————————————————————————————————————— ## +## ABUNDANCE OBSERVATION MODELS +## ———————————————————————————————————————————————————————————————————————————————————— ## +## QTYPE: +## 0 = FIXED SURVEY Q (specify log(mean) for prior log(mean)) +## 1 = CONSTANT Q (use MLE for q and optional informative prior) +## 2 = RANDOM WALK Q (use prior mean & sd for penalized random walk) +## ———————————————————————————————————————————————————————————————————————————————————— ## +1 # -number of surveys (n_it_nobs) + 2 # -QTYPE (see legend above) + 0 # -prior log(mean) + 0 # -prior sd (set to 0 for uniformative prior) + 1 # -Estimation phase +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + +## ———————————————————————————————————————————————————————————————————————————————————— ## +## OTHER MISCELANEOUS CONTROLS ## +## ———————————————————————————————————————————————————————————————————————————————————— ## 0 # 1 -verbose ADMB output (0=off, 1=on) 1 # 2 -recruitment model (1=beverton-holt, 2=ricker) 0.100 # 3 -std in observed catches in first phase. @@ -132,8 +204,13 @@ 0.50 # 13 -fraction of total mortality that takes place prior to spawning 0 # 14 -number of prospective years to add to syr. 0 # 15 -switch for IFD distribution in selectivity simulations -## -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## + + + + + +## ———————————————————————————————————————————————————————————————————————————————————— ## ## MARKER FOR END OF CONTROL FILE (eofc) -## ------------------------------------------------------------------------- ## +## ———————————————————————————————————————————————————————————————————————————————————— ## 999 \ No newline at end of file diff --git a/src/R/R4iSCAM.R b/src/R/R4iSCAM.R index 6908b9be..ee134c17 100644 --- a/src/R/R4iSCAM.R +++ b/src/R/R4iSCAM.R @@ -15,8 +15,8 @@ # | .RFILES <- List of R functions to source from the lib directory. # .PWD <- "/Users/stevenmartell1/Documents/iSCAM/examples/PacificHake14/R" .PWD <- "/Users/stevenmartell1/Documents/iSCAM-project/src/R" -.LIB <- "../../../dist/R/lib/" -.WIN <- "../../../dist/R/iScamWin2.txt" +.LIB <- "../../../src/R/lib/" +.WIN <- "../../../src/R/iScamWin2.txt" setwd(.PWD) .FIGUREDIR <- "../FIGS/" .RFILES <- list.files(.LIB,pattern="\\.[Rr]$") @@ -37,9 +37,6 @@ for(nm in .RFILES) source(file.path(.LIB, nm), echo=FALSE) .plotCatch( M ) .plotIndex( M ) .plotWeightAtAge( M ) -.plotAgeComps( M ) -.plotAgeCompResiduals( M ) -.plotAgeSummary( M ) .plotSpawnBiomass( M ) .plotDepletion( M ) .plotMortality( M ) @@ -47,9 +44,15 @@ for(nm in .RFILES) source(file.path(.LIB, nm), echo=FALSE) .plotStockRecruit( M ) .plotRecruitsPerSpawner( M ) .plotSurveyFit( M ) +.plotQ( M ) .plotCatchResidual( M ) .plotIndexResidual( M ) .plotRecruitmentResidual( M ) +.plotAgeComps( M ) +.plotAgeCompResiduals( M ) +.plotAgeSummary( M ) +.plotAgeBars( M ) +.plotSlx( M ) # |----------------------------------------------------------------------------------| diff --git a/src/R/lib/plotAgeBars.R b/src/R/lib/plotAgeBars.R new file mode 100644 index 00000000..2e676e19 --- /dev/null +++ b/src/R/lib/plotAgeBars.R @@ -0,0 +1,70 @@ +# plotAgeBars.R +library(reshape2) +library(ggplot2) + + +.plotAgeBars <- function( M ) +{ + n <- length( M ) + cat(".plotAgeBars\n") + id <- grep("d3_A[1-9]",names(M[[1]])) + iu <- grep("A_hat[1-9]",names(M[[1]])) + hr <- c("Year","Gear","Area","Group","Sex","AgeErr") + mdf <- NULL + mcol<- -1:-7 + for(i in 1:n) + { + getDF <- function(x) + { + ix <- id[x] + jx <- iu[x] + + age <- seq(M[[i]]$n_A_sage[x],M[[i]]$n_A_nage[x]) + df <- data.frame(V="Observed",M[[i]][[ix]]) + df[,mcol] <- df[,mcol]/rowSums(df[,mcol],na.rm=TRUE) + # cat("\nSum of observed = ",rowSums(df[,-1:-6])) + colnames(df) = c("V",hr,paste(age)) + + dp <- data.frame(V="Predicted",M[[i]][[ix]][,1:6],M[[i]][[jx]]) + dp[,mcol] <- dp[,mcol]/rowSums(dp[,mcol],na.rm=TRUE) + # cat("\nSum of predicted = ",rowSums(dp[,-1:-6])) + colnames(dp) = c("V",hr,paste(age)) + + dfp <- rbind(df,dp) + dfp <- subset(dfp,Year %in% seq(M[[i]]$syr,M[[i]]$nyr)) + return(dfp) + } + + B <- lapply(1:length(id),getDF) + } + + + + barplot <- function(B) + { + mB <- melt(B,id.vars=c("V",hr)) + mB$BroodYear <- mB$Year - as.integer(mB$variable) + + + O <- subset(mB,V=="Observed") + I <- subset(mB,V=="Predicted") + + klbl <- .GEAR[O$Gear] + albl <- .AREA[O$Area] + glbl <- .GROUP[O$Group] + slbl <- .SEX[O$Sex+1] + + p <- ggplot(O,aes(variable,as.double(value),fill=BroodYear)) + p <- p + geom_bar(stat="identity",alpha=0.7) + p <- p + scale_fill_distiller(type = "div") + p <- p + geom_line(data=I,aes(as.numeric(variable),value),color="red") + p <- p + labs(x="Age",y="Proportion") + p <- p + ggtitle(paste(klbl,":Area",albl,":Group",glbl,":Sex",slbl)) + p <- p + guides(title.position = "top") + p <- p + facet_wrap(~Year)+coord_flip() + return (p + .THEME + theme(legend.position="top") + theme_bw(8)) + } + + P <- lapply(B,barplot) + print(P) +} \ No newline at end of file diff --git a/src/R/lib/plotAgeComps.R b/src/R/lib/plotAgeComps.R index f6dbf33b..6545e601 100644 --- a/src/R/lib/plotAgeComps.R +++ b/src/R/lib/plotAgeComps.R @@ -1,64 +1,108 @@ # Rscript for plotting age-comp bubble plots. # Steven Martell # Aug 28, 2012 -require(ggplot2) -require(reshape2) +library(ggplot2) +library(reshape2) .plotAgeComps <- function( M ) { n <- length(M) cat(".plotAgeComps\n") + id <- grep("d3_A[1-9]",names(M[[1]])) mdf <- NULL for( i in 1:n ) { - A <- data.frame(M[[i]]$d3_A) - # Ensure proportions are being plotted. - A[,-1:-6] <- A[,-1:-6]/rowSums(A[,-1:-6],na.rm=TRUE) - age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) - # year gear area group sex - A <- data.frame(Model=names(M)[i],A) - colnames(A) <- c("Model","Year","Gear","Area","Group","Sex","AgeErr",paste(age)) - mdf <- rbind(mdf,A) + + getDF <- function(x) + { + ix <- id[x] + + df <- data.frame(M[[i]][ix]) + df <- data.frame(Model=names(M)[i],df) + age <- seq(M[[i]]$n_A_sage[x],M[[i]]$n_A_nage[x]) + colnames(df) <- c("Model","Year","Gear","Area","Group","Sex","AgeErr",paste(age)) + + return(df) + } + + B <- lapply(1:length(id),getDF) + + + + # A <- data.frame(M[[i]]$d3_A) + # # Ensure proportions are being plotted. + # A[,-1:-6] <- A[,-1:-6]/rowSums(A[,-1:-6],na.rm=TRUE) + # age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) + # # year gear area group sex + # A <- data.frame(Model=names(M)[i],A) + # colnames(A) <- c("Model","Year","Gear","Area","Group","Sex","AgeErr",paste(age)) + # mdf <- rbind(mdf,A) } - mdf <- melt(mdf,id.vars=c("Model","Year","Gear","Area","Group","Sex","AgeErr")) - BroodYear <- mdf$Year-as.double(mdf$variable) - mdf <- cbind(mdf,BroodYear) - print(head(mdf,3)) + mB <- melt(B,id.vars=c("Model","Year","Gear","Area","Group","Sex","AgeErr")) + BroodYear <- mB$Year-as.double(mB$variable) + mB <- cbind(mB,BroodYear) + + # mdf <- melt(mdf,id.vars=c("Model","Year","Gear","Area","Group","Sex","AgeErr")) + # BroodYear <- mdf$Year-as.double(mdf$variable) + # mdf <- cbind(mdf,BroodYear) + # print(head(mdf,3)) - p <- ggplot(mdf,aes(factor(Year),variable,size=value)) + p <- ggplot(mB,aes((Year),variable,size=value)) p <- p + geom_point(alpha=0.75,aes(colour=factor(BroodYear))) - p <- p + scale_size_area(max_size=10) + p <- p + scale_size_area(max_size=5) p <- p + labs(x="Year",y="Age",size="Count") - p <- p + facet_wrap(~Model+Sex+Gear,scales="free") + p <- p + facet_wrap(~L1+Model+Gear+AgeErr+Sex,scales="free") p <- p + scale_colour_discrete(guide="none") - print(p + .THEME) + print(p + .THEME + theme(legend.position="top")) } .plotAgeCompResiduals <- function( M ) { n <- length(M) cat(".plotAgeCompResiduals\n") + id <- grep("d3_A[1-9]",names(M[[1]])) + iu <- grep("A_nu[1-9]",names(M[[1]])) mdf <- NULL for( i in 1:n ) { - A <- cbind(M[[i]]$d3_A[,1:6],M[[i]]$A_nu) - A <- data.frame(A) - age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) - # year gear area group sex - A <- data.frame(Model=names(M)[i],A) - colnames(A) <- c("Model","Year","Gear","Area","Group","Sex","Err",paste(age)) - mdf <- rbind(mdf,subset(A,A$Year>=M[[i]]$syr)) + getDF <- function(x) + { + ix <- id[x] + jx <- iu[x] + + df <- data.frame(M[[i]][[ix]][,1:6],M[[i]][jx]) + df <- data.frame(Model=names(M)[i],df) + age <- seq(M[[i]]$n_A_sage[x],M[[i]]$n_A_nage[x]) + colnames(df) <- c("Model","Year","Gear","Area","Group","Sex","AgeErr",paste(age)) + + return(df) + } + + B <- lapply(1:length(id),getDF) + + + + # A <- cbind(M[[i]]$d3_A[,1:6],M[[i]]$A_nu) + # A <- data.frame(A) + # age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) + # # year gear area group sex + # A <- data.frame(Model=names(M)[i],A) + # colnames(A) <- c("Model","Year","Gear","Area","Group","Sex","Err",paste(age)) + # mdf <- rbind(mdf,subset(A,A$Year>=M[[i]]$syr)) } - mdf <- melt(mdf,id.vars=c("Model","Year","Gear","Area","Group","Sex","Err")) - print(head(mdf,3)) - p <- ggplot(mdf,aes(factor(Year),variable,col=factor(sign(value)),size=abs(value))) + mB <- melt(B,id.vars=c("Model","Year","Gear","Area","Group","Sex","AgeErr")) + + # mdf <- melt(mdf,id.vars=c("Model","Year","Gear","Area","Group","Sex","Err")) + # print(head(mdf,3)) + + p <- ggplot(mB,aes((Year),variable,col=factor(sign(value)),size=abs(value))) p <- p + geom_point(alpha=0.75) - p <- p + scale_size_area(max_size=6) + p <- p + scale_size_area(max_size=5) p <- p + labs(x="Year",y="Age",size="Residual",colour="Sign") - p <- p + facet_wrap(~Model+Sex+Gear,scales="free") - print(p + .THEME) + p <- p + facet_wrap(~L1+Model+Gear+AgeErr+Sex,scales="free") + print(p + .THEME + theme(legend.position="top")) } @@ -66,40 +110,77 @@ require(reshape2) { n <- length(M) cat(".plotAgeSummary\n") + id <- grep("d3_A[1-9]",names(M[[1]])) + iu <- grep("A_hat[1-9]",names(M[[1]])) mdf <- NULL + cat("length of M ",length(M),"\n") for( i in 1:n ) { - age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) - - # Predicted data - A <- cbind(M[[i]]$d3_A[,1:6],M[[i]]$A_hat) - A <- data.frame(A) - A[,-1:-6] <- A[,-1:-6]/rowSums(A[,-1:-6],na.rm=TRUE) - agA <- aggregate(A[,-1:-6],by=list(A[,2],A[,3],A[,4],A[,5]),FUN=mean,na.rm=TRUE) - colnames(agA) = c("Gear","Area","Group","Sex",paste(age)) - - # Observed data - O <- data.frame(M[[i]]$d3_A) - O[,-1:-6] <- O[,-1:-6]/rowSums(O[,-1:-6],na.rm=TRUE) - agO <- aggregate(O[,-1:-6],by=list(O[,2],O[,3],A[,4],O[,5]),FUN=mean,na.rm=TRUE) - colnames(agO) = c("Gear","Area","Group","Sex",paste(age)) - - # Create data frame - df <- rbind(cbind(type="Predicted",agA),cbind(type="Observed",agO)) - # year gear area group sex - df <- data.frame(Model=names(M)[i],df) - colnames(df) <- c("Model","Type","Gear","Area","Group","Sex",paste(age)) - mdf <- rbind(mdf,df) + + getDF <- function(x) + { + ix <- id[x] + jx <- iu[x] + + P <- data.frame(M[[i]][[ix]][,1:6],M[[i]][jx]) + P[,-1:-6] <- P[,-1:-6]/rowSums(P[,-1:-6],na.rm=TRUE) + age <- seq(M[[i]]$n_A_sage[x],M[[i]]$n_A_nage[x]) + aP <- aggregate(P[,-1:-6],by=list(P[,2],P[,3],P[,4],P[,5],P[,6]),FUN=mean,na.rm=TRUE) + colnames(aP) <- c("Gear","Area","Group","Sex","AgeErr",paste(age)) + + O <- data.frame(M[[i]][[ix]]) + O[,-1:-6] <- O[,-1:-6]/rowSums(O[,-1:-6],na.rm=TRUE) + age <- seq(M[[i]]$n_A_sage[x],M[[i]]$n_A_nage[x]) + aO <- aggregate(O[,-1:-6],by=list(O[,2],O[,3],O[,4],O[,5],O[,6]),FUN=mean,na.rm=TRUE) + colnames(aO) <- c("Gear","Area","Group","Sex","AgeErr",paste(age)) + + # create data frame + df <- rbind(data.frame(Type="Predicted",aP),data.frame(Type="Observed",aO)) + df <- data.frame(Model=names(M)[i],df) + colnames(df) <- c("Model","Type","Gear","Area","Group","Sex","AgeErr",paste(age)) + + return(df) + } + + B <- lapply(1:length(id),getDF) + + + + + # age <- seq(min(M[[i]]$n_A_sage),max(M[[i]]$n_A_nage)) + + # # Predicted data + # A <- cbind(M[[i]]$d3_A[,1:6],M[[i]]$A_hat) + # A <- data.frame(A) + # A[,-1:-6] <- A[,-1:-6]/rowSums(A[,-1:-6],na.rm=TRUE) + # agA <- aggregate(A[,-1:-6],by=list(A[,2],A[,3],A[,4],A[,5],A[,6]),FUN=mean,na.rm=TRUE) + # colnames(agA) = c("Gear","Area","Group","Sex","Err",paste(age)) + + # # Observed data + # O <- data.frame(M[[i]]$d3_A) + # O[,-1:-6] <- O[,-1:-6]/rowSums(O[,-1:-6],na.rm=TRUE) + # agO <- aggregate(O[,-1:-6],by=list(O[,2],O[,3],A[,4],O[,5],O[,6]),FUN=mean,na.rm=TRUE) + # colnames(agO) = c("Gear","Area","Group","Sex","Err",paste(age)) + + # # Create data frame + # df <- rbind(cbind(type="Predicted",agA),cbind(type="Observed",agO)) + # # year gear area group sex + # df <- data.frame(Model=names(M)[i],df) + # colnames(df) <- c("Model","Type","Gear","Area","Group","Sex","Err",paste(age)) + # mdf <- rbind(mdf,df) } - mdf <- melt(mdf,id.vars=c("Model","Type","Gear","Area","Group","Sex")) - print(head(mdf,3)) - p <- ggplot( mdf,aes(variable,value,col=Type) ) + mB <- melt(B,id.vars=c("Model","Type","Gear","Area","Group","Sex","AgeErr"),measured.vars=age) + + # mdf <- melt(mdf,id.vars=c("Model","Type","Gear","Area","Group","Sex","Err"),measured.vars=age) + # print(head(mdf,3)) + + p <- ggplot( mB,aes(variable,value,col=Type,shape=factor(Sex)) ) p <- p + geom_point(alpha=0.75) # p <- p + scale_area(range = c(0,10)) - p <- p + labs(x="Age",y="Mean proportion",colour="Type") - p <- p + facet_wrap(~Model+Sex+Gear,scales="free") - print(p + .THEME) + p <- p + labs(x="Age",y="Mean proportion",colour="Type",shape="Sex") + p <- p + facet_wrap(~L1+Model+Gear+AgeErr+Sex,scales="free") + print(p + .THEME + theme(legend.position="top")) } # .plotAgecomps <- function(repObj, meanAge = FALSE ) diff --git a/src/R/lib/plotCatchResiduals.R b/src/R/lib/plotCatchResiduals.R index fa260b30..80a48ed5 100644 --- a/src/R/lib/plotCatchResiduals.R +++ b/src/R/lib/plotCatchResiduals.R @@ -19,7 +19,7 @@ require(reshape2) } print(head(mdf,3)) - p <- ggplot(mdf,aes(x=factor(Year),Residual,fill=factor(Gear))) + p <- ggplot(mdf,aes(x=(Year),Residual,fill=factor(Gear))) p <- p + geom_bar(width=0.75,position="dodge",stat='identity') p <- p + labs(x="Year",y="log residual") p <- p + facet_wrap(~Model,scales="free") diff --git a/src/R/lib/plotQ.R b/src/R/lib/plotQ.R new file mode 100644 index 00000000..560241a4 --- /dev/null +++ b/src/R/lib/plotQ.R @@ -0,0 +1,41 @@ +# plotQ.R +# Steve Martell +# April 10, 2015 + + +.plotQ <- function ( M ) +{ + n <- length(M) + cat(".plotQ\n") + mdf <- NULL + for( i in 1:n ) + { + fit = M[[i]]$fit + yr = M[[i]]$d3_survey_data[,1] + kk = M[[i]]$d3_survey_data[,3] + + qt = M[[i]]$qt + qt = na.omit(as.vector(t(qt))) + ln.qt.mu = fit$est[grep("log_q_devs",fit$names)] + ln.qt.sd = fit$std[grep("log_q_devs",fit$names)] + ub = ln.qt.mu + 1.96 * ln.qt.sd + lb = ln.qt.mu - 1.96 * ln.qt.sd + + + df <- data.frame(Model=names(M)[i],year=yr,gear=kk,qt=qt) + names(df) <- c("Model","Year","Gear","Q") + df <- subset(df,Q>0) + + df <- cbind(df,"delta"=ln.qt.mu,"ub"=ub,"lb"=lb) + + mdf <- rbind(mdf,subset(df,Q>0)); + + } + + print(head(mdf)) + p <- ggplot(mdf) + geom_point(aes(Year,Q)) + p <- p + geom_point(aes(Year,delta)) + geom_pointrange(aes(Year,delta,ymax=ub,ymin=lb)) + p <- p + facet_grid(Gear~Model,scales="free") + print(p + .THEME) + +} \ No newline at end of file diff --git a/src/R/lib/plotRecruitmentResiduals.R b/src/R/lib/plotRecruitmentResiduals.R index ee5d53f7..5b139684 100644 --- a/src/R/lib/plotRecruitmentResiduals.R +++ b/src/R/lib/plotRecruitmentResiduals.R @@ -30,7 +30,7 @@ require(reshape2) print(head(mdf,3)) print(tail(mdf,3)) - p <- ggplot(mdf,aes(x=factor(Year),Residual)) + p <- ggplot(mdf,aes(x=(Year),Residual)) p <- p + geom_bar(width=0.75,position="dodge",stat='identity') p <- p + labs(x="Year",y="Recruitment (log residual)") p <- p + facet_wrap(~Model,scales="free") diff --git a/src/R/lib/plotSelex.R b/src/R/lib/plotSelex.R new file mode 100644 index 00000000..fa21f536 --- /dev/null +++ b/src/R/lib/plotSelex.R @@ -0,0 +1,44 @@ +# plotSelex.R +library(reshape2) +library(ggplot2) +library(dplyr) +library(gridExtra) + +.plotSlx <- function( M ) +{ + n <- length(M) + cat("plotSlx\n") + mdf <- NULL + decade <- seq(1900,2050,by=10) + for(i in 1:n) + { + df <- data.frame(Model=names(M)[i],logSel=M[[i]]$log_sel) + colnames(df) = c("Model","Gear","Area","Group","Sex","Year",M[[i]]$age) + df$Decade <- decade[findInterval(df$Year,decade)] + mdf <- rbind(mdf,melt(df,id=c("Model","Gear","Area","Group","Sex","Year","Decade"))) + } + + f <- unique(mdf$Area) + g <- unique(mdf$Group) + h <- unique(mdf$Sex) + k <- unique(mdf$Gear) + d <- unique(mdf$Decade) + ags <- as.matrix(expand.grid(k,f,g)) + + fnp <- function(ff) + { + dfx <- mdf %>% filter(Gear==ff[1]&Area==ff[2]&Group==ff[3]) + + p <- ggplot(dfx,aes(as.numeric(variable),exp(value),col=factor(Year))) + p <- p + geom_line() + p <- p + labs(x="Age",y="Relative Selectivity") + p <- p + ggtitle(paste("Gear",.GEAR[ff[1]])) + p <- p + facet_wrap(~Decade+Sex,scales="free_y") + return(p + .THEME + theme_bw(8) + guides(col=FALSE)) + } + + P <- apply(ags,1,fnp) + print(P) + + +} diff --git a/src/admb-code/include/compositionLikelihoods.hpp b/src/admb-code/include/compositionLikelihoods.hpp new file mode 100644 index 00000000..3e915ab2 --- /dev/null +++ b/src/admb-code/include/compositionLikelihoods.hpp @@ -0,0 +1,503 @@ +// negLogLikelihood.hpp +#include + +#ifndef negLogLikelihood_H +#define negLogLikelihood_H + +/** + * @file negLogLikelihood.hpp + * @defgroup Likelihoods + * @author Steven Martell + * @namespace acl + * @date Feb 10, 2014 + * @title Selectivity functions + * @details Uses abstract base class for computing negative loglikelihoods + * TODO: Still need to fix templates so it works with df1b2variables. + */ +namespace acl +{ + + /** + * @brief Base class for composition likelihoods. + * @details Virtual methods for negative loglikelihood and standardized residuals. + * + * @tparam DATA Data type argument. + * @tparam DVAR Variable type argument. + */ + template + class negLogLikelihood + { + private: + int r1,r2; /// index for first and last row. + int c1,c2; /// index for rist and last column. + ivector m_jmin; /// index for ragged start columns. + ivector m_jmax; /// index for ragged end columns. + imatrix m_jagg; /// index for aggregation among bins. + + DATA m_O; /// observed composition object + DATA m_rO; /// observed ragged composition object + + DVAR m_P; /// predicted composition object + DVAR m_rP; /// predicted ragged composition object + + + public: + // constructors + negLogLikelihood(){} + negLogLikelihood(const DATA& _O, const DVAR& _P) + :m_O(_O), m_P(_P) + { + r1 = m_O.rowmin(); + r2 = m_O.rowmax(); + c1 = m_O.colmin(); + c2 = m_O.colmax(); + + getRaggedVectors(); + } + + // pure virtual methods + virtual ~negLogLikelihood(){}; + virtual const dvariable nloglike() const = 0; + virtual const DVAR residual() const = 0; + + // virtual methods + inline + const DATA tailCompression(const DATA &_M) const; + + // methods + inline + void getRaggedVectors(double pmin = 0.0); + + inline + void aggregate(const double pmin = 0.0); + + inline + DVAR disaggregate(const DVAR &nu); + + // getters & setters + void set_O(DATA & O) { this -> m_O = O; } + DATA get_O() const{ return m_O; } + + DVAR get_P() const {return m_P; } + void set_P(DVAR &P) {this->m_P=P;} + + DATA get_rO() const {return m_rO; } + void set_rO(DATA &X){this->m_rO=X;} + + DVAR get_rP() const {return m_rP; } + void set_rP(DVAR &X){this->m_rP=X;} + + imatrix get_jagg() const {return m_jagg;} + + template + inline + const T compress(const T& _M) const + { + T R; + T M = _M; + R.allocate(r1,r2,m_jmin,m_jmax); + R.initialize(); + // fill ragged array R + for(int i = r1; i <= r2; i++ ) + { + M(i) /= sum(M(i)); + R(i)(m_jmin(i),m_jmax(i)) = M(i)(m_jmin(i),m_jmax(i)); + + // add cumulative sum to tails. + R(i)(m_jmin(i)) = sum(M(i)(c1,m_jmin(i))); + R(i)(m_jmax(i)) = sum(M(i)(m_jmax(i),c2)); + } + + return R; + } + + }; + + template + inline + void acl::negLogLikelihood::aggregate(const double pmin) + { + m_jmin.allocate(r1,r2); + m_jmax.allocate(r1,r2); + ivector n(r1,r2); + n.initialize(); + DATA oo(r1,r2,c1,c2); + DVAR pp(r1,r2,c1,c2); + + // get number of observations > pmin each year. + for(int i = r1; i <= r2; i++ ) + { + oo(i) = m_O(i)/sum(m_O(i)); + pp(i) = m_P(i)/sum(m_P(i)); + for(int j = c1; j <= c2; j++ ) + { + if( oo(i)(j) > pmin ) n(i)++; + } + } + + m_rO.allocate(r1,r2,1,n); + m_rP.allocate(r1,r2,1,n); + m_jagg.allocate(r1,r2,1,n); + m_rO.initialize(); + m_rP.initialize(); + m_jagg.initialize(); + + for(int i = r1; i <= r2; i++ ) + { + int k = 1; + for(int j = c1; j <= c2; j++ ) + { + m_rO(i)(k) += m_O(i)(j); + m_rP(i)(k) += m_P(i)(j); + if(k<=n(i)) m_jagg(i,k) = j; + if( oo(i)(j)>pmin && k + inline + DVAR acl::negLogLikelihood::disaggregate(const DVAR &nu) + { + // residuals to return. + DVAR NU; + NU.allocate(m_P); + NU.initialize(); + imatrix idx = this->get_jagg(); + for(int i = m_O.rowmin(); i <= m_O.rowmax(); i++ ) + { + for(int j = 1; j <= size_count(nu(i)); j++ ) + { + int jj = idx(i,j); /* code */ + NU(i)(jj) = nu(i)(j); + } + } + return NU; + } + + /** + * @brief Get column indicies for each row to compress matrix into ragged array. + * @details Determine the start and end positions of each row in which to compute + * the likelihoods where values < pmin are pooled into the tails of the distribution. + * + * @param pmin value to assume for the minium proportion. + * @return NULL Modifies private member variables: m_jmin and m_jmax. + */ + template + inline + void acl::negLogLikelihood::getRaggedVectors(double pmin) + { + m_jmin.allocate(r1,r2); + m_jmax.allocate(r1,r2); + + for(int i = r1; i <= r2; i++ ) + { + dvector o = m_O(i); + + m_jmin(i) = c1; // index for min column + m_jmax(i) = c2; // index for max column + dvector cumsum = o/sum(o); + for(int j = c1; j < c2; j++ ) + { + // min bin and plus bin column index. + cumsum(j) <= pmin ? m_jmin(i)++ : NULL; + 1.0 - cumsum(j) <= pmin ? m_jmax(i)-- : NULL; + cumsum(j+1) += cumsum(j); + + // o(j) > pmin ? n(i)++ : NULL; + } + } + + } + + + + + + + + + + + + + + + + + // |--------------------------------------------------------------------------------| + // | MULTIVARIATE LOGISTIC NEGATIVE LOGLIKELIHOOD derived class | + // |--------------------------------------------------------------------------------| + + /** + * @brief multivariate logistic negative log likelihood + * @details The negative log likeihood for the multivariate + * logistic distribution based on the MLE of the variance. + * + * @param NU Matrix of residuals + * @return the negative loglikelihood + */ + template + inline + const T dmvlogistic(const DVAR& NU) + { + int n = size_count(NU) - (NU.rowmax()-NU.rowmin()+1); + T var = 1.0/n * norm2(NU); + return n * log(var); + } + + template + inline + const DVAR dmvlogisticResidual(const DATA& O, const DVAR& P) + { + DVAR nu; + nu.allocate(P); + int i; + for( i = O.rowmin(); i<= O.rowmax(); i++) + { + nu(i) = log(O(i)/sum(O(i))) - log(P(i)/sum(P(i))); + nu(i) -= mean(nu(i)); + } + return nu; + } + + /** + * @brief Negative loglikelihood for compostion data using multivaritae logistic + * @author Steve Martell + * @param T usually a dmatrix + */ + template + class multivariteLogistic: public negLogLikelihood + { + private: + + DVAR m_nu; /// logistic residuals based on ragged object. + DVAR m_NU; /// residuals for rectangular matrix + + public: + // constructor + multivariteLogistic(const DATA &_O, const DVAR &_P, const double eps=0) + :negLogLikelihood(_O,_P) + { + // tail compression + DATA tmp = this->compress(this->get_O()) + eps; + this->set_rO(tmp); + DVAR vmp = this->compress(this->get_P()) + eps; + this->set_rP(vmp); + + + // residuals + this->aggregate(eps); + DVAR tnu = acl::dmvlogisticResidual(this->get_rO(),this->get_rP()); + set_nu(tnu); + + // residuals to return. + m_NU = this-> disaggregate(this -> get_nu()); + // m_NU.allocate(_P); + // m_NU.initialize(); + // imatrix idx = this->get_jagg(); + // for(int i = _O.rowmin(); i <= _O.rowmax(); i++ ) + // { + // for(int j = 1; j <= size_count(m_nu(i)); j++ ) + // { + // int jj = idx(i,j); /* code */ + // m_NU(i)(jj) = m_nu(i)(j); + // } + // } + } + + + const T nloglike() const + { + return acl::dmvlogistic(this->get_nu()); + } + + const DVAR residual() const + { + // This should return the residuals for the original container. + // Not the residuals for the ragged object. + + return this->get_NU(); + } + + DVAR get_nu() const {return m_nu; } + DVAR get_NU() const {return m_NU; } + void set_nu(DVAR &R){this->m_nu=R;} + + + }; + + + + + + + + + + + + + + + + + + + // |-------------------------------------------------------------------------------| + // | MULTINOMIAL DISTRIBUTION | + // |-------------------------------------------------------------------------------| + + /** + * @brief negative log density for multinomial distribution. + * @details Returns the negative log density for the multinomial distribution. + * Note that the sample size is based on the actual numbers in the matrix X. + * + * @param x Observed numbers in each bin + * @param p Predicted proportions in each bin. + * @param eps A tiny number to prevent log(0) + * @tparam T Templated variable (dvariable or df1b2variable) + * @return negative log density. + */ + template + T dmultinom(const DATA& x, const DVAR& p) + { + if(x.rowmin() != p.rowmin() || x.colmax() != p.colmax()) + { + cerr << "Index bounds do not macth in" + " dmultinom(const dvector& x, const dvar_vector& p)\n"; + exit(1); + } + + T ell = 0; + for(int i = x.rowmin(); i <= x.rowmax(); i++ ) + { + double n=sum(x(i)); + ell += -gammln(n+1.)+sum(gammln(x(i)+1.))-x(i)*log( p(i)/sum(p(i)) ); + } + return ell; + } + + + template + T dmultinom(const DATA& x, const DVAR& p, const T& log_vn) + { + if(x.rowmin() != p.rowmin() || x.colmax() != p.colmax()) + { + cerr << "Index bounds do not macth in" + " dmultinom(const dvector& x, const dvar_vector& p)\n"; + exit(1); + } + + T ell = 0; + T vn = exp(log_vn); + DVAR sobs = x; + for(int i = x.rowmin(); i <= x.rowmax(); i++ ) + { + + ell -= gammln(vn); + sobs(i) = vn * x(i) / sum(sobs(i)); + for(int j = x(i).indexmin(); j <= x(i).indexmax(); j++ ) + { + if( sobs(i,j) > 0.0 ) + ell += gammln(sobs(i,j)); + } + ell -= sobs(i) * log( p(i)/sum(p(i)) ); + } + return ell; + } + + template + inline + const DVAR dmultinomialResidual(const DATA& O, const DVAR& P) + { + // dvar_vector t1 = elem_div(o1-p1,sqrt(elem_prod(p1,1.-p1)/Nsamp)); + DVAR nu; + DVAR var; + nu.allocate(P); nu.initialize(); + var.allocate(P); var.initialize(); + + int i; + for( i = O.rowmin(); i<= O.rowmax(); i++) + { + nu(i) = O(i)/sum(O(i)) - P(i)/sum(P(i)); + var(i) = elem_prod(P(i),1.0-P(i)) / sum(O(i)); + nu(i) = elem_div(nu(i),sqrt(var(i))); + } + return nu; + } + + + template + class multinomial: public negLogLikelihood + { + private: + T m_nll; + DVAR m_nu; /// logistic residuals based on ragged object. + DVAR m_NU; /// residuals for rectangular matrix + + public: + // constructor with fixed sample size in _O + multinomial(const DATA &_O, const DVAR &_P, const double eps=0) + :negLogLikelihood(_O,_P) + { + // residuals + this -> aggregate(eps); + DVAR tnu = acl::dmultinomialResidual(this->get_rO(),this->get_rP()); + set_nu(tnu); + + // compute negative loglikelihood + T negLL = acl::dmultinom(this->get_rO(),this->get_rP(),eps); + set_nll(negLL); + + // residuals to return. + m_NU = this-> disaggregate(this -> get_nu()); + + } + + // constructor with estimated sample size in _O + multinomial(const DATA &_O, const DVAR &_P, const T &log_vn, const double eps=0) + :negLogLikelihood(_O,_P) + { + // residuals + this -> aggregate(eps); + DVAR tnu = acl::dmultinomialResidual(this->get_rO(),this->get_rP()); + set_nu(tnu); + + // compute negative loglikelihood + T negLL = acl::dmultinom(this->get_rO(),this->get_rP(),log_vn); + set_nll(negLL); + + + // residuals to return. + m_NU = this-> disaggregate(this -> get_nu()); + } + + const T nloglike() const + { + //return acl::dmultinom(this->get_rO(),this->get_rP()); + return get_nll(); + } + + const DVAR residual() const + { + return get_NU(); + } + + DVAR get_nu() const {return m_nu; } + DVAR get_NU() const {return m_NU; } + T get_nll() const {return m_nll; } + void set_nu(DVAR &R){this->m_nu=R;} + void set_nll(T &x) {this->m_nll=x; } + }; + +} + + + + + + + +#endif \ No newline at end of file diff --git a/src/admb-code/include/lib_iscam.h b/src/admb-code/include/lib_iscam.h index 650b961a..600579d5 100644 --- a/src/admb-code/include/lib_iscam.h +++ b/src/admb-code/include/lib_iscam.h @@ -8,6 +8,7 @@ #include "growthModels.hpp" #include "multinomial.h" #include "selex.hpp" +#include "compositionLikelihoods.hpp" #endif \ No newline at end of file diff --git a/src/admb-code/include/selex.hpp b/src/admb-code/include/selex.hpp index 7b3561c2..b50af4da 100644 --- a/src/admb-code/include/selex.hpp +++ b/src/admb-code/include/selex.hpp @@ -4,6 +4,8 @@ #define SELEX_HPP #include +#include +#include /** * @defgroup Selectivities @@ -21,6 +23,230 @@ */ namespace slx { + #define Register(x) slxInterface::template Init(); + + /** + * A Generic Selectivity Interface class. + */ + template + class slxInterface + { + private: + + typedef REAL_T(*EVAL_FUNCTION_PTR)(void* const); + + EVAL_FUNCTION_PTR eval_ptr; + + template static + REAL_T Evaluate_T(void* const pObj) + { + return static_cast (pObj)->Evaluate(); + } + + protected: + + // set a pointer to the evaluate function. + template void Init() + { + eval_ptr = (EVAL_FUNCTION_PTR) & Evaluate_T; + } + + public: + + slxInterface() + : eval_ptr(0){ + + } + + virtual ~ slxInterface(){ + + } + + inline REAL_T Evaluate(){ + assert(eval_ptr); // ensure Init() was called. + return (*eval_ptr)(this); + } + + }; + + + + + + + + /** + * Logistic Selectivity Class that inherits from the slxInterface + */ + template + class slx_Logistic: public slxInterface + { + private: + friend class slxInterface; + + // Evaluate should return a dvar_vector or a df1b2_vector + inline REAL_T Evaluate() { + //cout<<"Evaluating slx_Logistic"< + slx_Logistic(T1 &x, T2 &log_mu, T2 &log_sd) + :m_log_mu(log_mu),m_log_sd(log_sd),m_x(x) + { + Register(slx_Logistic); + } + }; + + + + + + + + /** + * Selectivity coefficients. + */ + + template + class slx_Coefficients: public slxInterface + { + private: + friend class slxInterface; + + inline REAL_T Evaluate() { + // cout<<"Evaluating slx_Coefficients"< + slx_Coefficients(T1 &x, T2 &log_sel_coeffs) + :m_x(x),m_log_sel_coeffs(log_sel_coeffs) + { + Register(slx_Coefficients); + } + }; + + + + + + + /** + * Cubic spline + */ + template + class slx_CubicSpline: public slxInterface + { + private: + friend class slxInterface; + + inline REAL_T Evaluate() { + //cout<<"Evaluating CubicSpline"< + slx_CubicSpline(T1 &x, T2 &log_spline_coeffs) + :m_x(x),m_log_spline_knots(log_spline_coeffs) + { + Register(slx_CubicSpline); + } + }; + + + + /** + * BiCubic spline + */ + template + class slx_BiCubicSpline: public slxInterface + { + private: + friend class slxInterface; + + inline REAL_T Evaluate() { + //cout<<"Evaluating BiCubicSpline"< + slx_BiCubicSpline(const T1 &x, const T1 &y, T2 &log_spline_coeffs, T2 &A) + :m_x(x),m_y(y),m_log_sel(A),m_log_spline_knots(log_spline_coeffs) + { + Register(slx_BiCubicSpline); + } + }; + + + + + + + + + + + + + template + const T plogis(const T &x, const T2 &mean, const T2 & sd) + { + return T2(1.0)/(T2(1.0)+exp(-(x-mean)/sd)); + } + /** * @ingroup Selectivities * @brief An abstract class for Selectivity functions. @@ -50,11 +276,6 @@ namespace slx { }; - template - const T plogis(const T &x, const T2 &mean, const T2 & sd) - { - return T2(1.0)/(T2(1.0)+exp(-(x-mean)/sd)); - } /** * @brief Logistic curve diff --git a/src/admb-code/iscam.tpl b/src/admb-code/iscam.tpl index ad4a550f..9efb0765 100644 --- a/src/admb-code/iscam.tpl +++ b/src/admb-code/iscam.tpl @@ -478,7 +478,7 @@ DATA_SECTION if(!mseFlag) { cout< Number of discrete selectivity blocks. // | + + /// April 13, issue 39. Changing the way selectivity/retention is controlled. + +!! #ifdef NEW_SELEX + init_imatrix slx_nBlocks(1,ngear,1,2); + int slx_nrow; + int ret_nrow + !! slx_nrow = sum(column(slx_nBlocks,1)); + !! ret_nrow = sum(column(slx_nBlocks,2)); + init_matrix slx_dControls(1,slx_nrow,1,13); + init_matrix ret_dControls(1,ret_nrow,1,13); + + ivector slx_nGearIndex(1,slx_nrow); /// index for fishing gear. + ivector slx_nIpar(1,slx_nrow); /// number of rows for each slx + ivector slx_nJpar(1,slx_nrow); /// number of cols for each slx + ivector slx_nSelType(1,slx_nrow); /// type of selectivity function + ivector slx_nAgeNodes(1,slx_nrow); /// number of age/size nodes + ivector slx_nYrNodes(1,slx_nrow); /// number of Year nodes + ivector slx_nSex(1,slx_nrow); /// index for sex (0=both, 1=female, 2=male) + ivector slx_phz(1,slx_nrow); /// phase of estimation or mirror index. + ivector slx_nsb(1,slx_nrow); /// start of block year. + ivector slx_neb(1,slx_nrow); /// end of block year. + vector slx_sel_mu(1,slx_nrow); + vector slx_sel_sd(1,slx_nrow); + vector slx_lam1(1,slx_nrow); + vector slx_lam2(1,slx_nrow); + vector slx_lam3(1,slx_nrow); + + + LOC_CALCS + slx_nGearIndex = ivector(column(slx_dControls,1)); + slx_nSelType = ivector(column(slx_dControls,2)); + slx_nSex = ivector(column(slx_dControls,5)); + slx_nAgeNodes = ivector(column(slx_dControls,6)); + slx_nYrNodes = ivector(column(slx_dControls,7)); + slx_phz = ivector(column(slx_dControls,8)); + slx_nsb = ivector(column(slx_dControls,12)); + slx_neb = ivector(column(slx_dControls,13)); + slx_sel_mu = column(slx_dControls,3); + slx_sel_sd = column(slx_dControls,4); + slx_lam1 = column(slx_dControls,9); + slx_lam2 = column(slx_dControls,10); + slx_lam3 = column(slx_dControls,11); + + // • Count number of selectivity parameters required for each slx_type + for(i = 1; i <= slx_nrow; i++) + { + slx_nIpar(i) = 1; + + switch(slx_nSelType(i)) + { + // • logistic selectivity + case 1: + slx_nJpar(i) = 2; + break; + + // • age-specific coefficients + case 2: + slx_nJpar(i) = int(nage - sage); + break; + + // • cubic spline over age/size + case 3: + slx_nJpar(i) = slx_nAgeNodes(i); + break; + + // • cubic spline over age/size each year + case 4: + slx_nJpar(i) = slx_nAgeNodes(i); + slx_nIpar(i) = slx_neb(i) - slx_nsb(i) + 1; + break; + + // • bicubic spline over age-size / years + case 5: + slx_nIpar(i) = slx_nAgeNodes(i); + slx_nJpar(i) = slx_nYrNodes(i); + break; + + // • logistic based on weight-at-age deviations. + case 6: + slx_nJpar(i) = 2; + break; + + // • logistic based on weight-at-age devs and scalar parameter. + case 7: + slx_nJpar(i) = 3; + break; + + // • age-specific coefficients between lb_age <= age <= ub_age + // case 8: + + // • logistic based on mean length-at-age. + case 11: + slx_nJpar(i) = 2; + break; + + // • size-specific coefficients. + //case 12: + + // • size-based cubic spline over mean size-at-age + case 13: + slx_nJpar(i) = slx_nAgeNodes(i); + break; + } + } + COUT(slx_nrow); + cout<<"Ok after new selex stuff in data section"< Standard deviation for constraint. // | m_nNodes -> number of nodes for cubic spline. // | m_nodeyear -> position of the nodes. - int nMdev; + int nMdev; init_int m_type; init_int Mdev_phz; - init_number m_stdev; + init_number m_stdev; init_int m_nNodes; init_ivector m_nodeyear(1,m_nNodes); @@ -1233,6 +1349,7 @@ DATA_SECTION init_ivector q_prior(1,nits); init_vector mu_log_q(1,nits); init_vector sd_log_q(1,nits); + init_ivector q_phz(1,nits); // |---------------------------------------------------------------------------------| // | Miscellaneous controls | @@ -1242,7 +1359,7 @@ DATA_SECTION // | 3 -> std in catch first phase // | 4 -> std in catch in last phase // | 5 -> assumed unfished in first year (0=FALSE, 1=TRUE) - // | 6 -> Maternal effects power parameter. + // | 6 -> Maternal effects power parameter (1=no maternal effects). // | 7 -> mean fishing mortality rate to regularize the solution // | 8 -> standard deviation of mean F penalty in first phases // | 9 -> standard deviation of mean F penalty in last phase. @@ -1291,12 +1408,13 @@ DATA_SECTION // | ilvec[2] -> number of surveys (nItNobs) // | ilvec[3] -> number of age-compisition data sets (nAgears) // | ilvec[4] -> container for recruitment deviations. - ivector ilvec(1,7); + ivector ilvec(1,8); !! ilvec = ngear; !! ilvec(1) = 1; !! ilvec(2) = nItNobs; - !! ilvec(3) = nAgears; - !! ilvec(4) = ngroup; + !! ilvec(3) = nItNobs; + !! ilvec(4) = nAgears; + !! ilvec(5) = ngroup; // |---------------------------------------------------------------------------------| @@ -1324,6 +1442,14 @@ DATA_SECTION if( dCatchData(i)(1) > nyr ) ft_count --; } + for( k = 1; k <= nItNobs; k++ ) + { + for( i = 1; i <= n_it_nobs(k); i++ ) + { + if( d3_survey_data(k)(i)(1) > nyr) qdev_count(k) --; + } + } + // Retrospective counter for n_A_nobs n_naa.initialize(); for( k = 1; k <= nAgears; k++ ) @@ -1343,6 +1469,7 @@ DATA_SECTION // | - adjust sel_blocks to new syr // | - Reduce ft_count so as not to bias estimates of ft. // | - Establish prospective counter for Composition data n_saa; + // | - Reduce qdev_count ivector n_saa(1,nAgears); @@ -1374,6 +1501,14 @@ DATA_SECTION if( iyr < syr ) n_saa(k)++; } } + + for( k = 1; k <= nItNobs; k++ ) + { + for( i = 1; i <= n_it_nobs(k); i++ ) + { + if( d3_survey_data(k)(i)(1) < syr) qdev_count(k) --; + } + } END_CALCS @@ -1408,29 +1543,16 @@ DATA_SECTION // | MANAGEMENT STRATEGY EVALUATION INPUTS // |---------------------------------------------------------------------------------| // | - - //LOC_CALCS - // ifstream ifile("Halibut2012.mse"); - // if(ifile) - // { - // cout<<"Vader is happy"< residuals between estimated R and R from S-R curve (process err) // | matrix log_rt(1,n_ag,syr-nage+sage,nyr); - matrix nlvec(1,7,1,ilvec); + matrix nlvec(1,8,1,ilvec); matrix epsilon(1,nItNobs,1,n_it_nobs); matrix xi(1,nItNobs,1,n_it_nobs); matrix it_hat(1,nItNobs,1,n_it_nobs); @@ -1659,7 +1814,7 @@ PARAMETER_SECTION matrix rt(1,ngroup,syr+sage,nyr); matrix delta(1,ngroup,syr+sage,nyr); - + // |---------------------------------------------------------------------------------| // | THREE DIMENSIONAL ARRAYS // |---------------------------------------------------------------------------------| @@ -1712,6 +1867,7 @@ PARAMETER_SECTION sdreport_vector sd_depletion(1,ngroup); sdreport_matrix sd_log_sbt(1,ngroup,syr,nyr+1); + PRELIMINARY_CALCS_SECTION @@ -1721,7 +1877,7 @@ PRELIMINARY_CALCS_SECTION // | - nf is a function evaluation counter. // | - SimFlag comes from the -sim command line argument to simulate fake data. // | - + nf=0; if( testMSY ) { @@ -1767,9 +1923,15 @@ RUNTIME_SECTION PROCEDURE_SECTION initParameters(); - + + #ifndef NEW_SELEX calcSelectivities(isel_type); - + #endif + + #ifdef NEW_SELEX + calcSelex(); + #endif + calcTotalMortality(); calcNumbersAtAge(); @@ -1833,6 +1995,9 @@ FUNCTION void calcSdreportVariables() sd_log_sbt(g) = log(sbt(g)); } + + + if( verbose ) { cout<<"**** Ok after calcSdreportVariables ****"< slx_nsb(k)?syr:slx_nsb(k); + int yr2 = nyr < slx_neb(k)?nyr:slx_neb(k); + int nn = slx_nIpar(k)-1; + + slx::slxInterface *ptrSlx[nn]; + for( i = 0; i <= nn; i++ ) + { + ptrSlx[i] = NULL; + } + slx::slxInterface *ptrSlxM = NULL; + + switch(slx_nSelType(k)) + { + // logistic selectivity based on age. + case 1: + for( j = 0; j < slx_nIpar(k); j++ ) + { + p1 = slx_log_par(k,j+1,1); + p2 = slx_log_par(k,j+1,2); + ptrSlx[j] = new slx::slx_Logistic(age,p1,p2); + } + break; + + // age-specific selectivity coefficients. + case 2: + for( j = 0; j < slx_nIpar(k); j++ ) + { + dvar_vector slx_theta = slx_log_par(k)(j+1); + ptrSlx[j] = new slx::slx_Coefficients(age,slx_theta); + } + break; + + // cubic spline + case 3: + for( j = 0; j < slx_nIpar(k); j++ ) + { + dvar_vector slx_theta = slx_log_par(k)(j+1); + ptrSlx[j] = new slx::slx_CubicSpline(age,slx_theta); + } + break; + + // • cubic spline over age/size each year + case 4: + for( j = 0; j < slx_nIpar(k); j++ ) + { + dvar_vector slx_theta = slx_log_par(k)(j+1); + ptrSlx[j] = new slx::slx_CubicSpline(age,slx_theta); + } + break; + + // • bicubic spline over age and year knots + case 5: + dvar_matrix tmp(yr1,yr2,sage,nage); + dvector iyr(1,slx_nYrNodes(k)); + dvector iag(1,slx_nAgeNodes(k)); + iyr.fill_seqadd(0,1.0/(slx_nYrNodes(k)-1)); + iag.fill_seqadd(0,1.0/(slx_nAgeNodes(k)-1)); + dvar_matrix slx_theta = slx_log_par(k); + tmp.initialize(); + + ptrSlxM = new slx::slx_BiCubicSpline(iag,iyr,slx_theta,tmp); + break; + } + + // fill arrays with selectivity coefficients. + // NOTES: + // • h = index for sex (0=both, 1=female, 2=male) + // • If slx_nSex(k) == 0, then apply same slx curve to both sexes. + // Do this by looping over area and group, and assign to specific sex. + j = 0; + int f,g,h; + int h_sex = slx_nSex(k); + int kgear = slx_nGearIndex(k); + for(int ig = 1; ig <= n_ags; ig++ ) + { + f = n_area(ig); + g = n_group(ig); + h = n_sex(ig); + + // if h_sex == 0, then you need to reset j = 0 + if ( h_sex == 0 ) j = 0; + + // if !h_sex, skip the process if current group is not the right sex + if ( h_sex != 0 && h != h_sex) continue; + int igrp = pntr_ags(f,g,h); + // Fill vectors of selex + if (ptrSlx[j]) + { + for(i = yr1; i <= yr2; i++) + { + log_sel(kgear)(igrp)(i) = ptrSlx[j] -> Evaluate(); + if(slx_nSelType(k) == 4 && j < slx_nIpar(k)) j++; + } + } + + // Fill matrix of selex + if (ptrSlxM) + { + log_sel(kgear)(igrp).sub(yr1,yr2) = ptrSlxM -> Evaluate(); + } + + //subtract mean to ensure mean(exp(log_sel))==1 + for(i = yr1; i <= yr2; i++) + { + log_sel(kgear)(igrp)(i) -= log( mean(mfexp(log_sel(kgear)(ig)(i))) ); + } + } + + if( !ptrSlxM ) delete ptrSlxM; + if( !*ptrSlx ) delete *ptrSlx; + + + } + + if(verbose==1) cout<<"End of CalcSelex"< nyr) continue; + // total numbers-at-age + // ta.initialize(); + // for( h = 1; h <= nsex; h++ ) + // { + // ig = pntr_ags(f,g,h); + // ta += N(ig)(i); + // } + if( h ) // age comps are sexed (h > 0) { ig = pntr_ags(f,g,h); @@ -2579,10 +2875,12 @@ FUNCTION calcComposition { dvar_vector pred_ca = ca * age_age(e); A_hat(kk)(ii) = pred_ca(n_A_sage(kk),n_A_nage(kk)); + + // predicted plus group age. if( n_A_nage(kk) < nage ) - { - A_hat(kk)(ii)(n_A_nage(kk)) += sum( pred_ca(n_A_nage(kk)+1,nage) ); - } + { + A_hat(kk)(ii)(n_A_nage(kk)) += sum( pred_ca(n_A_nage(kk)+1,nage) ); + } } else { @@ -2761,6 +3059,7 @@ FUNCTION calcTotalCatch fishery that would remove potential spawn that would not be surveyed. - d3_survey_data: (iyr index(it) gear area group sex wt timing) - for MLE of survey q, using weighted mean of zt to calculate q. + - March 30, 2015. Added deviation in q for random walk. TODO list: [] - add capability to accomodate priors for survey q's. @@ -2773,11 +3072,12 @@ FUNCTION calcSurveyObservations int ii,kk,ig,nz; double di; - dvariable ftmp; + dvariable ftmp,zbar; dvar_vector Na(sage,nage); dvar_vector va(sage,nage); dvar_vector sa(sage,nage); epsilon.initialize(); + qt.initialize(); xi.initialize(); it_hat.initialize(); @@ -2839,42 +3139,90 @@ FUNCTION calcSurveyObservations dvar_vector t1 = rowsum(V); - dvar_vector zt = log(it) - log(t1(iz,nz)); - dvariable zbar = sum(elem_prod(zt,wt)); - q(kk) = mfexp(zbar); - - // | survey residuals - epsilon(kk).sub(iz,nz) = elem_div(zt - zbar,it_log_se(kk)(iz,nz)); - it_hat(kk).sub(iz,nz) = q(kk) * t1(iz,nz); + dvar_vector zt = log(it); + - // | SPECIAL CASE: penalized random walk in q process error only. - if( q_prior(kk)==2 ) + // | March 30, 2015. Issue #37. + // | Added switch for 3 q_prior(kk) options. + // | q_prior(kk) = 0 = constant fixed Q + // | q_prior(kk) = 1 = constant MLE Q + // | q_prior(kk) = 2 = penalized random walk in Q + switch( q_prior(kk) ) { - epsilon(kk).initialize(); - dvar_vector fd_zt = first_difference(zt); - dvariable zw_bar = sum(elem_prod(fd_zt,wt(iz,nz-1))); - epsilon(kk).sub(iz,nz-1) = elem_div(fd_zt - zw_bar,it_log_se(kk)(iz,nz-1)); - qt(kk)(iz) = exp(zt(iz)); - for(ii=iz+1;ii<=nz;ii++) - { - qt(kk)(ii) = qt(kk)(ii-1) * exp(fd_zt(ii-1)); - } - it_hat(kk).sub(iz,nz) = elem_prod(qt(kk)(iz,nz),t1(iz,nz)); + case 0: // Constant fixed Q + q(kk) = exp(mu_log_q(kk)); + it_hat(kk).sub(iz,nz) = q(kk) * t1(iz,nz); + zt -= log(it_hat(kk).sub(iz,nz)); + break; + + case 1: // Constant MLE Q + zt -= log(t1(iz,nz)); + zbar = sum(elem_prod(zt,wt)); + q(kk) = mfexp(zbar); + + // | survey residuals + it_hat(kk).sub(iz,nz) = q(kk) * t1(iz,nz); + zt -= zbar; + //epsilon(kk).sub(iz,nz) = elem_div(zt,it_log_se(kk)(iz,nz)); + break; + + case 2: // Penalized random walk in Q + int jj = 1; // index for qdevs. + zt -= log(t1(iz,nz)); + qt(kk)(iz) = exp( zt(iz) + log_q_devs(kk)(jj) ); + for(ii=iz+1; ii<=nz; ii++) + { + jj ++; + qt(kk)(ii) = qt(kk)(ii-1) * exp(log_q_devs(kk)(jj)); + + } + it_hat(kk).sub(iz,nz) = elem_prod(qt(kk)(iz,nz),t1(iz,nz)); + zt -= log(qt(kk)(iz,nz)); + + //epsilon(kk).sub(iz,nz)= elem_div(zt,it_log_se(kk)(iz,nz)); + break; } - // | MIXED ERROR MODEL for random walk in q - if( q_prior(kk)==3 ) + // Standardized observation error residuals. + epsilon(kk).sub(iz,nz) = elem_div(zt,it_log_se(kk)(iz,nz)); + + // Standardized process error residuals. + if(active(log_q_devs(kk))) { - dvar_vector proerr = zt - zbar; - qt(kk)(ii) = exp(zbar + proerr(iz)); - for(ii=iz+1;ii<=nz;ii++) - { - proerr(ii) = zt(ii) - zt(ii-1); - qt(kk)(ii) = qt(kk)(ii-1) * exp(proerr(ii)); - } - xi(kk).sub(iz,nz) = elem_div(proerr,it_log_pe(kk)(iz,nz)); - it_hat(kk).sub(iz,nz) = elem_prod(qt(kk)(iz,nz),t1(iz,nz)); + //dvar_vector fd_qt = first_difference( log_q_devs(kk) ); + dvar_vector fd_qt = first_difference( log(qt(kk)(iz,nz)) ); + xi(kk).sub(iz,nz-1) = elem_div(fd_qt.shift(iz),it_log_pe(kk)(iz,nz-1)); } + +// TO BE DEPRECATED +// // | SPECIAL CASE: penalized random walk in q process error only. +// if( q_prior(kk)==2 ) +// { +// epsilon(kk).initialize(); +// dvar_vector fd_zt = first_difference(zt); +// dvariable zw_bar = sum(elem_prod(fd_zt,wt(iz,nz-1))); +// epsilon(kk).sub(iz,nz-1) = elem_div(fd_zt - zw_bar,it_log_se(kk)(iz,nz-1)); +// qt(kk)(iz) = exp(zt(iz)); +// for(ii=iz+1;ii<=nz;ii++) +// { +// qt(kk)(ii) = qt(kk)(ii-1) * exp(fd_zt(ii-1)); +// } +// it_hat(kk).sub(iz,nz) = elem_prod(qt(kk)(iz,nz),t1(iz,nz)); +// } +// +// // | MIXED ERROR MODEL for random walk in q +// if( q_prior(kk)==3 ) +// { +// dvar_vector proerr = zt - zbar; +// qt(kk)(ii) = exp(zbar + proerr(iz)); +// for(ii=iz+1;ii<=nz;ii++) +// { +// proerr(ii) = zt(ii) - zt(ii-1); +// qt(kk)(ii) = qt(kk)(ii-1) * exp(proerr(ii)); +// } +// xi(kk).sub(iz,nz) = elem_div(proerr,it_log_pe(kk)(iz,nz)); +// it_hat(kk).sub(iz,nz) = elem_prod(qt(kk)(iz,nz),t1(iz,nz)); +// } } if(verbose)cout<<"**** Ok after calcSurveyObservations ****"< multivariate logistic using conditional MLE of the variance for weight. // | - 2 -> multnomial, assumes input sample size as n in n log(p) - // | - 3 -> logistic normal w no autocorrelation. + // | - 3 -> logistic normal w no autocorrelation. // | - Both likelihoods pool pmin (d_iscamCntrl(16)) into adjacent yearclass. // | - PSEUDOCODE: // | - => first determine appropriate dimensions for each of nAgears arrays (naa) @@ -3120,6 +3469,15 @@ FUNCTION calcObjectiveFunction // | TODO: // | [ ] - change A_nu to data-type variable, does not need to be differentiable. // | [ ] - issue 29. Fix submatrix O, P for prospective analysis & sex/area/group. + + // Testing new abstract compositionLikelihood class. SUXS. + + //acl::negLogLikelihood *ptr_AgeCompLike[nAgears-1]; + + + + + A_nu.initialize(); for(k=1;k<=nAgears;k++) { @@ -3145,9 +3503,11 @@ FUNCTION calcObjectiveFunction P(ii) = A_hat(k)(i).sub(n_A_sage(k),n_A_nage(k)); ii ++; } - //if( iyr <= nyr ) naa++; - //if( iyr < syr ) iaa++; } + + + + //dmatrix O = trans(trans(d3_A_obs(k)).sub(n_A_sage(k),n_A_nage(k))).sub(iaa,naa); //dvar_matrix P = trans(trans(A_hat(k)).sub(n_A_sage(k),n_A_nage(k))).sub(iaa,naa); @@ -3159,20 +3519,42 @@ FUNCTION calcObjectiveFunction logistic_student_t cLST_Age( O,P,dMinP(k),dEps(k) ); switch( int(nCompLikelihood(k)) ) { - case 1: - nlvec(3,k) = dmvlogistic(O,P,nu,age_tau2(k),dMinP(k)); + case 1: // multivariate Logistic + nlvec(4,k) = dmvlogistic(O,P,nu,age_tau2(k),dMinP(k)); + // cout<<"like: "<(O,P,dMinP(k)); + // nlvec(4,k) = ptr_AgeCompLike[k-1] -> nloglike(); + // nu = ptr_AgeCompLike[k-1] -> residual(); + break; - case 2: - nlvec(3,k) = dmultinom(O,P,nu,age_tau2(k),dMinP(k)); + + case 2: // multinomial with fixed sample size. + //ptr_AgeCompLike[k-1] = new acl::multinomial(O,P,dMinP(k)); + //nlvec(4,k) = ptr_AgeCompLike[k-1] -> nloglike(); + //nu = ptr_AgeCompLike[k-1] -> residual(); + //COUT(nlvec(4,k)); + nlvec(4,k) = dmultinom(O,P,nu,age_tau2(k),dMinP(k)); + //COUT(nlvec(4,k)); + //exit(1); break; + + case 6: // Multinomial with estimated effective sample size. + //ptr_AgeCompLike[k-1] = new acl::multinomial(O,P,log_degrees_of_freedom(k),dMinP(k)); + //nlvec(4,k) = ptr_AgeCompLike[k-1] -> nloglike(); + //nu = ptr_AgeCompLike[k-1] -> residual(); + // deprecate + nlvec(4,k) = mult_likelihood(O,P,nu,log_degrees_of_freedom(k)); + + break; + case 3: if( !active(log_age_tau2(k)) ) // LN1 Model { - nlvec(3,k) = cLN_Age(); + nlvec(4,k) = cLN_Age(); } else { - nlvec(3,k) = cLN_Age( exp(log_age_tau2(k)) ); + nlvec(4,k) = cLN_Age( exp(log_age_tau2(k)) ); } // Residual @@ -3182,16 +3564,17 @@ FUNCTION calcObjectiveFunction age_tau2(k) = cLN_Age.get_sigma2(); } break; + case 4: //logistic_normal cLN_Age( O,P,dMinP(k),dEps(k) ); if( active(phi1(k)) && !active(phi2(k)) ) // LN2 Model { - nlvec(3,k) = cLN_Age(exp(log_age_tau2(k)),phi1(k)); + nlvec(4,k) = cLN_Age(exp(log_age_tau2(k)),phi1(k)); } if( active(phi1(k)) && active(phi2(k)) ) // LN3 Model { - nlvec(3,k) = cLN_Age(exp(log_age_tau2(k)),phi1(k),phi2(k)); + nlvec(4,k) = cLN_Age(exp(log_age_tau2(k)),phi1(k),phi2(k)); } // Residual @@ -3206,11 +3589,11 @@ FUNCTION calcObjectiveFunction case 5: // Logistic-normal with student-t if( !active(log_degrees_of_freedom(k)) ) { - nlvec(3,k) = cLST_Age(); + nlvec(4,k) = cLST_Age(); } else { - nlvec(3,k) = cLST_Age(exp(log_degrees_of_freedom(k))); + nlvec(4,k) = cLST_Age(exp(log_degrees_of_freedom(k))); } // Residual @@ -3220,14 +3603,12 @@ FUNCTION calcObjectiveFunction age_tau2(k) = cLST_Age.get_sigma2(); } break; - case 6: // Multinomial with estimated effective sample size. - nlvec(3,k) = mult_likelihood(O,P,nu,log_degrees_of_freedom(k)); - break; case 7: // Multivariate-t - nlvec(3,k) = multivariate_t_likelihood(O,P,log_age_tau2(k), + nlvec(4,k) = multivariate_t_likelihood(O,P,log_age_tau2(k), log_degrees_of_freedom(k), phi1(k),nu); age_tau2(k) = exp(value(log_age_tau2(k))); + cout<<"Here we go dude"<slx_nsb(kr)?syr:slx_nsb(kr); + int yr2 = nyr 0) + { + nlvec(7,k) += slx_lam2(k) * square(diff); + } + } + } + } + } + + // penalty in time-varying changes (nlvec(8,k)). + if( slx_nSelType(kr) == 4 || slx_nSelType(kr) == 5) + { + for(ig=1;ig<=n_ags;ig++) + { + dvar_matrix trans_log_sel = trans( log_sel(k)(ig) ); + for(j=sage;j<=nage;j++) + { + dvar_vector df1 = first_difference(trans_log_sel(j)); + dvar_vector df2 = first_difference(df1); + nlvec(8,k) += slx_lam3(k)/(nage-sage+1)*norm2(df2); + } + } + } + } + + } + #endif + + #ifndef NEW_SELEX for(k=1;k<=ngear;k++) { - if(active(sel_par(k))) + + if( active(sel_par(k)) ) { //if not using logistic selectivity then //CHANGED from || to && May 18, 2011 Vivian @@ -3286,12 +3727,12 @@ FUNCTION calcObjectiveFunction { //curvature in selectivity parameters dvar_vector df2 = first_difference(first_difference(log_sel(k)(ig)(i))); - nlvec(5,k) += lambda_1(k)/(nage-sage+1)*df2*df2; + nlvec(6,k) += lambda_1(k)/(nage-sage+1)*df2*df2; //penalty for dome-shapeness for(j=sage;j<=nage-1;j++) if(log_sel(k,ig,i,j)>log_sel(k,ig,i,j+1)) - nlvec(6,k)+=lambda_2(k) + nlvec(7,k)+=lambda_2(k) *square( log_sel(k,ig,i,j)-log_sel(k,ig,i,j+1) ); } } @@ -3303,7 +3744,11 @@ FUNCTION calcObjectiveFunction Mar 13, 2013, added 2nd difference penalty on isel_type==5 */ - if( isel_type(k)==4 || isel_type(k)==5 || n_sel_blocks(k) > 1 ) + if( lambda_3(k) && + (isel_type(k)==4 || + isel_type(k)==5 || + n_sel_blocks(k) > 1) + ) { for(ig=1;ig<=n_ags;ig++) { @@ -3312,20 +3757,55 @@ FUNCTION calcObjectiveFunction for(j=sage;j<=nage;j++) { dvar_vector df2 = first_difference(first_difference(trans_log_sel(j))); - nlvec(7,k) += lambda_3(k)/(nage-sage+1)*norm2(df2); + nlvec(8,k) += lambda_3(k)/(nage-sage+1)*norm2(df2); } } } } } - + #endif // |---------------------------------------------------------------------------------| // | CONSTRAINTS FOR SELECTIVITY DEVIATION VECTORS // |---------------------------------------------------------------------------------| // | [?] - TODO for isel_type==2 ensure mean 0 as well. // | + + #ifdef NEW_SELEX + for(k=1;k<=slx_nrow;k++) + { + if( active(slx_log_par(k)) && slx_nSelType(k)!=1 ) + { + dvariable s = 0; + //bicubic spline version ensure col mean = 0 + if(slx_nSelType(k)==5) + { + dvar_matrix tmp = trans(slx_log_par(k)); + for(j=1;j<=tmp.rowmax();j++) + { + s=mean(tmp(j)); + lvec(1)+=10000.0*s*s; + } + } + + if( slx_nSelType(k)==2 || + slx_nSelType(k)==3 || + slx_nSelType(k)==4 || + slx_nSelType(k)==12 ) + { + dvar_matrix tmp = slx_log_par(k); + for(j=1;j<=tmp.rowmax();j++) + { + s=mean(tmp(j)); + lvec(1)+=10000.0*s*s; + } + } + } + } + #endif + + #ifndef NEW_SELEX for(k=1;k<=ngear;k++) { if( active(sel_par(k)) && @@ -3358,6 +3838,7 @@ FUNCTION calcObjectiveFunction } } } + #endif // |---------------------------------------------------------------------------------| @@ -3413,7 +3894,7 @@ FUNCTION calcObjectiveFunction qvec.initialize(); for(k=1;k<=nits;k++) { - if(q_prior(k) == 1 ) + if(q_prior(k) == 1 && sd_log_q(k) != 0 ) { qvec(k) = dnorm( log(q(k)), mu_log_q(k), sd_log_q(k) ); } @@ -4605,11 +5086,32 @@ REPORT_SECTION { REPORT(n_A_sage); REPORT(n_A_nage); + + for(k = 1; k<=nAgears; k++) + { + adstring lbl = "d3_A"+str(k); + report<