# Computing SE for PSD indices

Demonstrate how to compute SE for PSD indices.
PSD
Size Structure
Author

Derek H. Ogle

Published

Dec 9, 2016

Modified

Dec 19, 2022

Note

The following packages are loaded for use below.

library(FSA)   ## for psdAdd(), rcumsum(), headtail()
library(dplyr) ## for mutate(), filter()

A reader of Ogle (2016) asked how to compute standard errors (SE) for the various PSD indices with the usual SE equation of a proportion (i.e., $$\sqrt{\frac{p(1-p)}{n}}$$. I demonstrate this calculation in this post

Warning

I caution against using these standard errors to produce confidence intervals for PSD-X values. If confidence intervals are the end product for this analysis, then I suggest using the binomial or multinomial distribution methods described in Ogle (2016).

I will use the same Inch Lake data from the Size Structure chapter of Ogle (2016), where I used psdAdd() (from FSA) to add a variable with the Gabelhouse length categories and eliminated all sub-stock size fish.1

1 this and next step is from this script) from Ogle (2016).

inchAlls <- read.csv("http://derekogle.com/IFAR/scripts/InchLake1113.csv") |>
filter(gcat!="substock") |>
droplevels()
#R|  No known Gabelhouse (PSD) lengths for: Iowa Darter
headtail(inchAlls)
#R|       gearType year       species  tl      gcat
#R|  1        fyke 2012 Black Crappie 246   quality
#R|  2        fyke 2012 Black Crappie 165     stock
#R|  3        fyke 2012 Black Crappie 330 memorable
#R|  1270     fyke 2013  Yellow Perch 348 memorable
#R|  1271     fyke 2013  Yellow Perch 135     stock
#R|  1272     fyke 2013  Yellow Perch 137     stock

I then computed the PSD-X values as shown in Ogle (2016). These will be the proportions (i.e., $$p$$) in the SE calculations.

( freq <- xtabs(~species+gcat,data=inchAlls) )
#R|                   gcat
#R|  species           stock quality preferred memorable
#R|    Black Crappie      33      40         8        81
#R|    Bluegill          571     199       223        43
#R|    Largemouth Bass    26      36         0         0
#R|    Pumpkinseed         5       2         0         0
#R|    Yellow Perch        4       0         0         1
iPSDs <- prop.table(freq,margin=1)*100
( PSDs <- t(apply(iPSDs,MARGIN=1,FUN=rcumsum)) )
#R|                   gcat
#R|  species           stock  quality preferred memorable
#R|    Black Crappie     100 79.62963  54.93827 50.000000
#R|    Bluegill          100 44.88417  25.67568  4.150579
#R|    Largemouth Bass   100 58.06452   0.00000  0.000000
#R|    Pumpkinseed       100 28.57143   0.00000  0.000000
#R|    Yellow Perch      100 20.00000  20.00000 20.000000

The row sums from the freq table are the total number of fish (that are stock size or greater) and will be $$n$$ in the SE calculations.

( sums <- rowSums(freq) )
#R|    Black Crappie        Bluegill Largemouth Bass     Pumpkinseed    Yellow Perch
#R|              162            1036              62               7               5

Finally, with these two results the SE for PSD-Q may be computed.

p <- PSDs[,"quality"]
SEs <- sqrt(p*(100-p)/sums)
round(SEs,1)
#R|    Black Crappie        Bluegill Largemouth Bass     Pumpkinseed    Yellow Perch
#R|              3.2             1.5             6.3            17.1            17.9

And can be repeat for other size indices; e.g., for PSD-P.

p <- PSDs[,"preferred"]
SEs <- sqrt(p*(100-p)/sums)
round(SEs,1)
#R|    Black Crappie        Bluegill Largemouth Bass     Pumpkinseed    Yellow Perch
#R|              3.9             1.4             0.0             0.0            17.9

Note

The SE formula is usually used for proportions (and then multiplied by 100 to get a SE for percentages). However, it can be shown algebraically to apply to percentages as well.

## References

Ogle, D. H. 2016. Introductory Fisheries Analyses with R. CRC Press, Boca Raton, FL.

## Citation

BibTeX citation:
@online{h. ogle2016,
author = {H. Ogle, Derek},
title = {Computing {SE} for {PSD} Indices},
date = {2016-12-09},
url = {https://fishr-core-team.github.io/fishR//blog/posts/2016-12-9_PSD_SE},
langid = {en}
}