@@ -3,14 +3,112 @@ gb_area <- ecodata::epu_sf |>
33 sf :: st_area() /
44 1e6 # convert to km2
55
6- ecodata :: ppr | >
6+ new_ppr <- ecodata :: ppr | >
77 dplyr :: filter(EPU == " GB" , Var == " Ryther" | Var == " PP" ) | >
88 dplyr :: select(- Units ) | >
99 tidyr :: pivot_wider(names_from = Var , values_from = Value ) | >
1010 dplyr :: mutate(
11- total_energy_required = (Ryther * gb_area ) | > units :: drop_units(),
12- ratio = total_energy_required / PP
13- ) | >
11+ total_fish = (Ryther * gb_area ) | > units :: drop_units(),
12+ ppr = (total_fish / 9 ) * (1 / 0.15 )^ 2.5 ,
13+ ratio = ppr / PP
14+ )
15+
16+ new_ppr | >
17+ ggplot2 :: ggplot(ggplot2 :: aes(x = Time , y = ratio )) +
18+ ggplot2 :: geom_point() +
19+ ggplot2 :: geom_line() +
20+ ggplot2 :: theme_bw()
21+
22+ # compare to pp from annual_chl_pp
23+
24+ ecodata :: plot_annual_chl_pp(
25+ report = " NewEngland" ,
26+ EPU = " GB" ,
27+ varName = " pp" ,
28+ plottype = " total"
29+ )
30+
31+ fix <- dplyr :: filter(
32+ dplyr :: mutate(
33+ ecodata :: annual_chl_pp ,
34+ Time = as.integer(gsub(" A_" , " " , Time ))
35+ ),
36+ EPU == " GB" ,
37+ Var == " PPD_ANNUAL_MTON"
38+ )
39+
40+ ppr <- ecodata :: ppr | >
41+ dplyr :: filter(EPU == " GB" , Var == " PP" )
42+
43+ combined_pp <- tibble :: tibble(
44+ Time = c(
45+ fix $ Time ,
46+ ppr | >
47+ dplyr :: pull(Time )
48+ ),
49+ pp = c(
50+ fix $ Value ,
51+ ppr | >
52+ dplyr :: pull(Value )
53+ ),
54+ source = c(
55+ rep(" annual_chl_pp" , nrow(fix )),
56+ rep(
57+ " ppr" ,
58+ length(
59+ ppr | >
60+ dplyr :: pull(Time )
61+ )
62+ )
63+ )
64+ )
65+
66+ combined_pp | >
67+ ggplot2 :: ggplot(ggplot2 :: aes(x = Time , y = pp , color = source )) +
68+ ggplot2 :: geom_point() +
69+ ggplot2 :: geom_line() +
70+ ggplot2 :: theme_bw()
71+
72+ # compare catch in ppr to ecodata::comdat
73+
74+ ecodata :: comdat | >
75+ dplyr :: filter(Var == " Landings" , EPU == " GB" ) | >
76+ ggplot2 :: ggplot(ggplot2 :: aes(x = Time , y = Value )) +
77+ ggplot2 :: geom_point() +
78+ ggplot2 :: geom_line() +
79+ ggplot2 :: theme_bw() +
80+ ggplot2 :: geom_point(
81+ data = new_ppr ,
82+ ggplot2 :: aes(y = total_fish ),
83+ color = " red"
84+ )
85+
86+ reconstructed_ppr <- dplyr :: full_join(
87+ ecodata :: comdat | >
88+ dplyr :: filter(Var == " Landings" , EPU == " GB" ) | >
89+ dplyr :: select(Time , Value ) | >
90+ dplyr :: mutate(Var = " Landings" ),
91+ ecodata :: annual_chl_pp | >
92+ dplyr :: filter(EPU == " GB" , Var == " PPD_ANNUAL_MTON" ) | >
93+ dplyr :: select(Time , Value ) | >
94+ dplyr :: mutate(
95+ Var = " PPD_ANNUAL_MTON" ,
96+ Time = as.integer(gsub(" A_" , " " , Time ))
97+ )
98+ ) | >
99+ tidyr :: pivot_wider(names_from = Var , values_from = Value ) | >
100+ dplyr :: mutate(
101+ ppr = (Landings / 9 ) * (1 / 0.15 )^ 2.5 ,
102+ # fill in missing pp
103+ PPD_ANNUAL_MTON = ifelse(
104+ is.na(PPD_ANNUAL_MTON ),
105+ mean(PPD_ANNUAL_MTON , na.rm = TRUE ),
106+ PPD_ANNUAL_MTON
107+ ),
108+ ratio = ppr / PPD_ANNUAL_MTON
109+ )
110+
111+ reconstructed_ppr | >
14112 ggplot2 :: ggplot(ggplot2 :: aes(x = Time , y = ratio )) +
15113 ggplot2 :: geom_point() +
16114 ggplot2 :: geom_line() +
0 commit comments