Fo es Ecology and Managemen 549 (2023) 121501
A ailable online 19 Oc obe 2023
0378-1127/© 2023 The Au ho s. Published by Else ie B.V. This is an open access a icle unde he CC BY-NC-ND license (h p://c ea i ecommons.o g/licenses/by-
nc-nd/4.0/).
Base-age in a ian models o p edic ing indi idual ee accumula ed
annual esin yield using wo apping me hods in ma i ime pine (Pinus
pinas e Ai .) o es s in no h-wes e n Spain
´
Osca L´
opez-´
Al a ez
*
, Luis F anco-V´
azquez, Manuel Ma ey-Pe ez
Resea ch G oup PROePLA GI-1716, Depa men o Plan P oduc ion and Enginee ing P ojec s, Highe Poly echnic School o Enginee ing, Campus Te a, Uni e si y o
San iago de Compos ela, 27002 Lugo, Spain
ARTICLE INFO
Keywo ds:
Pine esin
P edic ing
Algeb aic di e ence app oach
Dynamic equa ion
Be alan y-Richa ds
ABSTRACT
In sou he n Eu ope, especially in Spain and Po ugal, ma i ime pine esin is one o he main non- imbe o es
p oduc s. A e su e ing a c isis a he end o he 20 h cen u y, i is cu en ly a g owing sec o . In Spain,
depending on he a ea, he managemen o pine o es s is one o he pilla s o he na ional bioeconomy. In
addi ion o imbe p oduc ion, hese o es s may be o ien ed owa ds esin p oduc ion only, o esin p oduc ion
as a complemen a y ac i i y o imbe p oduc ion. In bo h cases, as in any sec o , i is essen ial o ha e ools o
manage and an icipa e p oduc ion, especially in he new con ex o he bioeconomy. Fo his eason, he aim o
his s udy is o de elop a dynamic model o es ima e he accumula ed esin yield du ing he esin p oduc ion
season. Fo his s udy, 180 ees om h ee plo s loca ed in he no hwes o he Ibe ian Peninsula we e esin
apped using wo ex ac ion me hods (non-mechanized and mechanized ci cula no ching) and s imulan pas es.
Fou base models we e used om which eigh equa ions we e de i ed using ADA and GADA echniques. The
mos e icien equa ions, bo h o modelling wi h he ain da a and o p edic ion wi h he es da a, we e hose
de i ed om he Be alan y-Richa ds model. The RRMSE was 23% o he non-mechanised me hod and 29% o
he mechanised ci cula me hod. The esul s o his s udy make i possible o add he cumula i e annual esin
yield o ma i ime pine o he p ocesses ha he Be alan y-Richa ds equa ion is capable o modelling.
Fu he mo e, he g ea e sa ili y o hese models will be o g ea use o he o es manage in op imising he
annual ha es ing season as well as o he scien i ic communi y.
1. In oduc ion
Th oughou his o y, he esou ces gene a ed by o es s ha e p o-
ided di e en ypes o goods ele an o socie y (Sheppa d e al., 2020).
The social and en i onmen al bene i s associa ed wi h mul i- unc ional
o es managemen p ac ices ha e an impo an impac , especially in
u al a eas (Lo i´
c e al., 2022). Wo ldwide, he main o es esou ce in
e ms o olume and e enue is imbe , bu o he commodi ies such as
non- imbe o es p oduc s (NTFPs) a e also aluable (Sa deshpande and
Shackle on, 2019). Acco ding o FAO (2015), he NTFPs a e “goods
de i ed om o es s ha a e angible and physical objec s o biological
o igin o he han wood” like ood addi i es, ib es, esins o gums. The
high added alue o NTFPs makes hem a pe ec complemen o imbe
p oduc ion, p o iding an ex a sou ce o income be ween imbe ha -
es s (He n´
andez-Rod íguez e al., 2017). One such bene i , de i ed
om hei na u e as a biop oduc , is hei ole in mi iga ing he ca bon
oo p in and educing he use o ossil-based ma e ials (Solomon,
2016). This is due o he abili y o some o hem ( ib es, gum, o esin) o
ix CO
2
and o eplace pe oleum de i a i es (Demko and Macha a,
2022).
In Eu ope, acco ding o Vacik e al. (2020), he economic impo ance
o plan s-de i ed NTFPs emains modes , wi h he epo ed alue o
NTFPs eaching 1.6 billion eu os. Da a o s a is ics a e di icul o
ob ain, since a la ge p opo ion o NTFPs a e in ended o sel -
consump ion (Lo i´
c e al., 2020; Winkel e al., 2022). Despi e his,
he e a e Eu opean coun ies ha dedica e exclusi e esou ces o his
ype o ac i i ies, such as F ance, Po ugal, and Spain, which a e
commi ed o NTFPs wi h he highes added alue (Wol slehne e al.,
2019). In Spain, he alue o NTFPs in ela ion o imbe a ies be ween
10 % and 50 %, depending on he sou ce (Díaz-Bal ei o e al., 2020).
* Co esponding au ho .
E-mail add esses: [email p o ec ed] (´
O. L´
opez-´
Al a ez), [email p o ec ed] (L. F anco-V´
azquez), [email p o ec ed] (M. Ma ey-Pe ez).
Con en s lis s a ailable a ScienceDi ec
Fo es Ecology and Managemen
jou nal homepage: www.else ie .com/loca e/ o eco
h ps://doi.o g/10.1016/j. o eco.2023.121501
Recei ed 19 June 2023; Recei ed in e ised o m 7 Sep embe 2023; Accep ed 13 Oc obe 2023
Fo es Ecology and Managemen 549 (2023) 121501
2
Acco ding o MITECO (2022), he main Spanish NTFPs p oduc s a e
u les and edible ungi (106.6 million eu os), co k (83.6 million eu os),
ches nu (18 million eu os), pine esin (10.3 million eu os) and pine nu s
(2.6 million eu os). Pine esin is he main de ensi e ba ie agains
ex e nal a acks caused by pes s and diseases o pine ees (Rissanen
e al., 2019). Spanish pine esin p oduc ion has expe ienced ups and
downs h oughou his o y. Since he 1980 s, esin apping ac i i y has
los impo ance un il i p ac ically ceased a he beginning o his cen-
u y (Soli˜
no e al., 2018). Bu since 2000 s i has been g owing in
impo ance again wi hin he na ional o es y sec o and is expec ed o
con inue o g ow in he coming yea s (G´
omez-Ga cía e al., 2022). This
is e lec ed in he inc ease in publica ions on he subjec in some egions
o Spain in ecen yea s, as in he case o Galicia (Ga cía-M´
eijome e al.,
2023; L´
opez-´
Al a ez e al., 2023; Touza e al., 2021; V´
azquez-Gonz´
alez
e al., 2021; Zas e al., 2020a,b). This g ow h is due o i s mul iple in-
dus ial uses, he en i onmen al bene i s associa ed wi h ex ac ion and
he boos he sec o is ecei ing om public and p i a e ini ia i es
(Soli˜
no e al., 2018; Oono e al., 2020).
The p og essi e decline o he sec o in he pas cen u y has led o a
lack o echnical expe ise and knowledge gene a ion, which does no
co espond o o he sec o s such as he wood sec o (G´
eno a e al., 2014;
L´
opez-´
Al a ez e al., 2023a). In addi ion, exis ing da abases ha e
p oblems in ob aining comple e da a se ies (Sainz e al, 2010; Calama
e al., 2020). The complexi y o collec ing hese da a se ies is due o he
la ge numbe o in e media e measu emen s ha need o be made and
he ac ha i is a e y labo ious p ocess. In o de o make he sec o
mo e compe i i e, he wo k o p oduce s mus be acili a ed by
p o iding hem wi h di e en ypes o esou ces. Rega ding issues
ela ed o he mechaniza ion o apping p ocess, a ious au ho s ha e
compa ed he pe o mance be ween di e en mechanized and non-
mechanized ex ac ion me hods (Rod íguez-Ga cía e al., 2016; L´
opez-
Al a ez e al., 2023b; Ga cía-M´
eijome e al. 2023) and be ween di e en
s imulan pas es used (Neis e al., 2018; Micha ila e al., 2021). In e ms
o planning ools, he e is a lack o s a is ically obus models ha can
p edic inal p oduc ion be o e o once he apping season has s a ed.
This ype o modelling is essen ial o mul i unc ional o es esou ce
planning, as i allows he s ock and expec ed alue o esinous pine o be
es ima ed.
This lack o obus models may be due o he complexi y o iden i-
ying he explana o y a iables o o al p oduc ion, as i depends on he
in e ac ion o mul iple in e ela ed a iables (Zas e al., 2020a). The
pape by L´
opez-Al a ez e al. (2023b) shows he impo an biome ic
and p oduc ion a iabili y o pines wi hin and be ween plo s, and p e-
ious s udies mainly ela e esin yield in empe a e zones o wo en i-
onmen al ac o s, empe a u e and wa e s ess (Zas e al., 2020b).
These esul s a e in line wi h he ma kedly seasonal beha iou o annual
esin p oduc ion, wi h he highes p oduc ions ob ained in he wa me
and less ainy mon hs (Rod íguez-Ga cía e al., 2015). This pa e n
causes he indi idual ee accumula ed esin p oduc ion cu e in he
season o ollow a sigmoidal cu e, a pa e n widely s udied in mul iple
ields (Cagla e al., 2018). Fo example, ee and s and g ow h beha es
simila ly when ep esen ed o e ime (Pio esan and Biondi, 2021).
Fo es e s use he cons uc ion o si e index cu es o model his end
o e ime. They use h ee main echniques: guide cu e, pa ame e
p edic ion and algeb aic di e ence app oach (Clu e e al., 1983). The
las me hodology has he ad an age ha i allows he de elopmen o
dynamic equa ions de i ed om de eloping age-in a ian models om
epea ed measu emen s aken o e ime (Jo dan e al., 2006). The main
me hods o de eloping hese equa ions a e he algeb aic di e ence
app oach (ADA), o malized by Bailey and Clu e (1974), and he
gene alized algeb aic di e ence app oach (GADA), in oduced by
Cieszewski and Bailey (2000). The ADA me hod has he limi a ion ha
he models de i ed om i only allow he c ea ion o amilies o
anamo phic cu es o hose wi h a single asymp o e (Cieszewski, 2001).
To sol e he limi a ions o he ADA me hod, he GADA me hodology was
de eloped, which allows he c ea ion o amilies o polymo phic cu es
wi h mul iple asymp o es (Cieszewski, 2002). In GADA, o ob ain hese
amilies o cu es, he base equa ions a e expanded acco ding o se e al
cha ac e is ics ha a e oo ed on di e en g ow h heo ies (Di´
eguez-
A anda e al., 2005). Acco ding o hese heo ies, he cu es should s a
om he o igin, be polymo phic, ha e mul iple asymp o es, and ha e
he same es ima ed alue be ween p edic ion and e e ence age (Cies-
zewski and Bailey, 2000).
The aim o his pape is o e alua e he easibili y o modelling he
indi idual ee annual accumula ed esin p oduc ion de eloping an
in a ian model o pine o es s in no hwes e n Spain. To ca y ou he
main objec i e, he pape has been s uc u ed in o wo pa s: (1) e al-
ua ion o he alidi y o he ADA and GADA me hods o modelling he
indi idual ee accumula ed annual esin yield and (2) he o ecas ing
pe o mance o models.
2. Ma e ial and me hods
2.1. S udy a ea and da a
The da a used o ca y ou he s udy we e limi ed due o he di icul y
o collec ion and he lack o p e ious s udies. As a esul , in o ma ion o
he da ase came om h ee newly es ablished plo s loca ed in low-
densi y pu e s ands o Pinus pinas e in he no hwes e n Spanish e-
gion o Galicia (Fig. 1). The s udy a ea is cha ac e ized by mode a e
annual a e age empe a u es (11–14 ◦C) and accumula ed annual
ain all be ween 800 and 1500 mm.
Du ing he mon hs o June o Oc obe 2021, a o al o 180 ees we e
selec ed om he plo s o esin apping p ocedu es. The minimum,
mean, and maximum diame e s a b eas heigh (dbh), o al heigh s (h),
and esin yield (Y) alues o he ees in each plo a e p esen ed in
Table 1. The a e age dbh mee s he legal equi emen s o esin p o-
cedu es in he egion. As a p elimina y s ep, he dasome ic a iables
we e es ed o hei lack o in luence on esin yield by calcula ing
Spea man co ela ions be ween esin and dbh (co =0.18, p- alue =
0.02) and h (co =0.07, p- alue =0.37), as nei he a iable ollowed a
no mal dis ibu ion acco ding o he Shapi o-Wilk es (dbh =W =0.96,
p- alue =0.001; h =W =0.98, p- alue =0.03).
I was also checked whe he he yields o he plo s belonged o he
same popula ion. Fo his pu pose, he Shapi o-Wilk es (W =0.969, p-
Fig. 1. Dis ibu ions a ea o he ma i ime pine s ands used in he s udy.
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
3
alue =0.001) was used o e i y ha he da a o yields as a unc ion o
plo did no sa is y he assump ion o no mali y, bu he Le ene es (F =
0.246, p- alue =0.781) did sa is y he assump ion o homoscedas ici y.
Ha ing e i ied his, he K uskal-Wallis es (Fig. 2) was used o es
whe he he e we e s a is ically signi ican di e ences be ween he
yields o each plo , which showed ha he e we e no s a is ically sig-
ni ican di e ences be ween hem.
Resin apping was ca ied ou using wo di e en me hods, applying
each me hod o 90 ees (30 ees pe plo using each me hod, a o al o
180 ees). One was a non-mechanized a ian o he “Ame ican”
ex ac ion me hod desc ibed in Rod íguez-Ga cía e al. (2014), used
mainly in he Ibe ian Peninsula and F ance. The o he ex ac ion me hod
was he mechanized ci cula g oo e me hod (L´
opez-´
Al a ez e al.,
2023b), in which a ci cula no ch is made in he s em o he ee and a
ci cula de ice is placed o collec he esin in a closed plas ic bag. The
non-mechanized me hod cu a 16 cm long and 3 cm wide ec angula
s ip o ba k and cambium, while he mechanized ci cula g oo e
me hod used a 5 cm diame e ool and a ba e y-powe ed sc ewd i e
(cu ing ci cum e ence o 15.71 cm). A new s ip was made e e y 14
days. By he end o he apping season, which las ed i e mon hs, a o al
o 10 cu s had been made. Two di e en pas es we e used as a esin
p oduc ion s imulan , e hephon (8 % e hephon (60 % / ), 14 % sul-
phu ic acid (50 % / ), 55 % dis illed wa e , 1.7 % polyso ba e, 1 %
ce yl alcohol, 4 % aseline, 5.5 % silica, 10.8 sawdus ) and ASACIF (1 %
salicylic acid, 25 % sulphu ic acid (96 % / ), 5 % p opylene glycol, 19
% whea s aw, 50 % dis illed wa e ). These pas es we e applied each
ime a new s ip was made. Due o incon eniences in he da a collec ion
p ocess, he e was a small g oup o ees o which some in e media e
yield could no be measu ed a he p ecise momen when he new cu -
ing was made and we e he e o e accumula ed in he nex weighing. To
i he models, pe iodic yields we e accumula ed mon hly, gi ing a o al
o i e measu emen s o each ee.
2.2. Model de elopmen
To model he beha iou o esin yield o e he yea , we chose o spli
he da a sample by ex ac ion me hod, bu no by s imulan pas e. This
da a spli ing was done because he e was a high a iabili y be ween he
Table 1
Dasome ic and esin yield cha ac e iza ion o he plo s u ilized in he s udy. dbh: diame e a b eas heigh (cm); h: o al heigh (m); Y: ee esin yield (g).
Plo dbh
min
dbh dbh
max
h
min
h h
max
Y
min
Y Y
max
Culle edo 21.4 31.6 45.9 12.0 16.1 21.0 407 1419.11 3359
Pan ´
on 23.1 32.6 47.9 17.0 22.5 27.0 642 1564.47 2862
Godos 25.5 39.5 54.4 17.3 20.4 24.1 339 1492.27 2974
Fig. 2. S a is ical di e ences es be ween he esin yields o each plo . The s a is ical es s be ween g oups we e pe o m wi h s a is ical so wa e R (R Co e Team,
2022) and he “ggs a splo ” package (Pa il, 2021).
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
4
yields in he da a se and s a is ically signi ican di e ences be ween he
ex ac ion me hods, as can be seen in Fig. 3. The e a e no s a is ical
di e ences be ween he s imulan pas es (Fig. 3).
2.2.1. Models es ed
Se e al models ha e been es ed, including bo h ADA and GADA
models de i ed om di e en base equa ions. The e o e, he e a e
equa ions wi h one o mo e a iable pa ame e s depending on he
speci ic si e. The gene al no a ion (Di´
eguez-A anda e al., 2005;
S´
anchez-Gonz´
alez e al., 2008; Panik, 2014) used o de elop he models
has been, a
1
, a
2
,…, a
n
o deno e he pa ame e s in he base models, and
b
1
, b
2
,…, b
n
o he pa ame e s ha we e i ed in he di e en o mu-
la ions o he base equa ions. Hence, he gene al o mula ion o he
de i ed models has he o m Y = (
0
,
1
, Y
0
, b
1
, b
2
, …, b
n
). These ypes
o equa ions ha e mul iple uses, such as desc ibing he g ow h o ees
and hei pheno ypic elemen s, o he g ow h o co k hickness, as
men ioned be o e. Al hough he e a e o he uses besides hose pu ely
o es y, like i s use in ishe ies (Flinn and Midway, 2021), mussels
(Fuen es-San os e al., 2017) o chicken (Ma a-Es ada e al., 2020)
g ow h models.
Di e en ADA and GADA o mula ions om base equa ions we e
es ed, bu inally he selec ed base equa ions ha achie ed con e gence
(Table 2) we e he Be alan y-Richa ds (Be alan y, 1949, 1957;
Richa ds, 1959), Ko (ci ed in Lundq is (1957)), Hoss eld (Hoss eld,
1822) and Weibull (Weibull, 1951; Yang e al., 1978). A o al o eigh
equa ions de i ed om he p e ious base equa ions we e i ed. Fou o
hem (Eqs. E1 o E4) ha e been based on he Be alan y-Richa ds
model. O which, h ee (Eqs. E1 o E3) we e de i ed by he ADA
me hodology, hus possessing only one si e-speci ic pa ame e . Equa ion
E1 gene a es amilies o anamo phic cu es wi h mul iple asymp o es,
whe eas Equa ions E2 and E3 gene a e polymo phic cu es wi h a single
asymp o e (Manso e al., 2021). The esul o de i ing he Be alan y-
Richa ds equa ion wi h he GADA me hodology is Eq. E4, i has mo e
han one pa ame e ha is si e-speci ic and allows gene a ing amilies o
polymo phic cu es wi h mul iple asymp o es (P ada e al., 2019). In he
case o Eqs E5 o E8, which a e de i ed om he a ious equa ions
men ioned abo e, hese a e models ob ained using he ADA me hodol-
ogy, only ha ing a si e-speci ic pa ame e .
2.2.2. Pa ame e es ima ion, model selec ion and alida ion
To pe o m he i and subsequen alida ion o he models, he da a
wi hou ou lie s was andomly s a i ied in o qua iles acco ding o he
ime lag p edic o a iable in aining (80 %) and es (20 %). The i s
se was used o adjus he model pa ame e s simul aneously (local and
global pa ame e s) using he dummy a iable me hod desc ibed in
Cieszewski e al. (2000), and selec he one wi h he bes goodness-o - i
s a is ics. The second se was used o e alua e he pe o mance o he
selec ed model in an unbiased manne .
Model selec ion was pe o med by nume ical and g aphical analysis
o he esiduals o he i s. To e alua e he pe o mance o he models on
he aining da a, goodness-o - i s a is ics such as.
pseudo-R
2
(R
2
), calcula ed acco ding o Schabenbe ge and Pie ce
(2001) using he esidual sum o squa es and he o al co ec ed sum o
squa es, Roo Mean Squa e E o (RMSE) and Akaike’s in o ma ion
c i e ion (AIC) we e used. The p edic i e pe o mance o he models was
e alua ed using he Rela i e Roo Mean Squa e E o (RRMSE). To
p o ide addi ional model de ails, we ha e included he bias ( he a e age
o he p edic ion esiduals) and MAE ( he a e age o he absolu e alues
o he p edic ion esiduals).
The s a is ical so wa e R (R Co e Team, 2022) was used o he
calcula ions. The “ sample” package (F ick e al., 2022) was used o
andomly sepa a e he sample in o es and aining se s. The “s a s”
package (R Co e Team, 2022) was used o pe o m he non-linea models
i ing. The ollowing unc ions and packages we e used o calcula e he
goodness-o - i s a is ics and measu e he accu acy o he models on he
es da a: he “aomisc” package (Ono i, 2020) o R
2
, he “Me ics”
package (Hamne and F asco, 2018) o RMSE, he “s a s” package (R
Co e Team, 2022) o AIC, and he “me ica” package (Co endo e al.,
2022) o RRMSE.
3. Resul s
3.1. Accumula ed esin yield pa e n
The smoo hed end o he esin yield ob ained a indi idual ha es s
du ing he collec ion o esin apping da a o each o he ex ac ion
me hods, ollows he pa e ns shown in Fig. 4. In he non-mechanized
Fig. 3. S a is ical di e ences be ween a) he ex ac ion me hods and b) and c) he s imulan pas es acco ding o he ex ac ion me hod. The s a is ical es s be ween
g oups we e pe o m wi h s a is ical so wa e R (R Co e Team, 2022) and he “ggs a splo ” package (Pa il, 2021).
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
5
apping me hod, which was he mos p oduc i e in e ms o o al esin
yield, he immedia e esin yield be ween s ips was inc eased om he
beginning o he season eaching he highes alues in he la e summe
mon hs and hen dec eased un il he end o he apping season. In
con as , he mechanised ci cula g oo e me hod showed a mo e
consis en pa e n. I s end inc eased om he s a o he esin
ex ac ion campaign un il 30 days a e he s a , emained s able un il
40 days be o e he end, and hen dec eased un il he end o he season.
Consequen ly, he ci cula me hod did no show as p onounced a peak as
he non-mechanised me hod.
As a consequence o he phenomenon shown in Fig. 4, when he
cumula i e yields o e he esin apping season o indi idual ees a e
plo ed as a unc ion o esin ex ac ion me hod, hey ollow a sigmoidal
pa e n (Fig. 5). I can also be seen ha he slope o he cu es o he non-
mechanised me hod is g ea e han he slope o he mechanised ci cula
me hod, which is a consequence o he same phenomenon.
3.2. Final models
The es ima ed pa ame e s, s anda d e o s o he es ima es,
goodness-o - i s a is ics and pe o mance on alida ion da a a e shown
in Table 3 and Table 4. All he i ed pa ame e s we e s a is ically sig-
ni ican a a 1 % le el.
3.3. Non-mechanized apping me hod
All he i ed models exhibi ed sa is ac o y pe o mance a e
adjus ing o cumula i e esin yield (Table 3). Upon e alua ion using he
aining da a, he R
2
alues ob ained anged om 0.818 o 0.873,
depending on he speci ic model. The RMSE alues a ied om 246.980
o 300.758, while he AIC alues anged om 17964.02 o 18474.63. As
a measu e o he good pe o mance o he models, we looked o hose
wi h he lowes RMSE and AIC alues. On he es da a, he models
demons a ed ai RRMSE alues below 30 % (Li e al., 2013; Despo o ic
e al., 2016), wi h none exceeding 27.018 % and he lowes being
23.189 %. The bias ange was be ween −24.34 and 19.74 g while he
MAE was be ween 164.38 and 177.78 g.
Du ing he aining phase, he op pe o me s in e ms o R
2
, RMSE
and AIC we e h ee models de i ed om he Be alan y-Richa ds
model, ollowed by one model de i ed om he Ko , Hoss eld, and
Weibull equa ions. Con e sely, he leas e ec i e models we e he
anamo phic and mul iple asymp o ic models de i ed om he
Be alan y-Richa ds equa ion, as well as he Ko model when sol ing
o he “a
1
″
pa ame e . Simila esul s we e ob ained when assessing he
p edic i e pe o mance o he models on he es da a.
Among he models, equa ion E3, which co esponds o he ADA
model when sol ing o he “a
3
″
pa ame e in he Be alan y-Richa ds
equa ion, exhibi ed he highes R
2
alue, as well as he lowes RMSE
and AIC alues. When assessing he p edic i e pe o mance o his
model on he es da a, an RRMSE o 23 %, a bias o −19.48 g and a MAE
o 170.38 g we e ob ained. Fu he mo e, g aphical analysis o he e-
siduals o he E3 model e ealed no e idence o he e oscedas ici y o
au oco ela ion be ween he esiduals (Fig. 6).
To ep esen he e olu ion o he cumula i e p oduc ions gene a ed
by he model (Fig. 6), a ime ( ) o 56 days om he beginning o he
apping campaign was selec ed. This pe iod was chosen as a e e ence
since i is wo mon hs a e he s a o he campaign and one mon h
be o e he hal way poin o he campaign, allowing o an ea ly es i-
ma ion o he inpu s o he season. To ep esen he p oduc ion cu es,
he di e ence be ween he maximum and minimum cumula i e p o-
duc ion in he ime was di ided by ou , and his amoun was added
Table 2
Base models and ADA and GADA o mula ions conside ed.
Base model Pa ame e ela ed o si e Solu ion o X wi h ini ial alues (
0
, Y
0
) Dynamic equa ion
Be alan y-Richa ds
Y=a1⋅(1−e−a2⋅ )a3
a
1
=X X0=Y0
(1−e−b2⋅ 0)b3 Y=Y0⋅(1−e−b2⋅
1−e−b2⋅ 0)b3 (E1)
a
2
=X
X0=
−ln(1− (Y0/b1)1/b3)
0 Y=b1⋅⎛
⎝1−(1−(Y0
b1)1/b3) / 0⎞
⎠
b3 (E2)
a
3
=X X0=ln(Y0/b1)
ln(a−e−b2⋅ 0)Y=b1⋅(Y0
b1)
ln(1−e−b2⋅ )
ln(1−e−b2⋅ 0)
(E3)
a
1
=e
x
a
3
=b
2
+b
3
/X X0=1
2⋅(lnY0−b2⋅L0+
(lnY0−b2⋅L0)2−4⋅b3⋅L0
√)
wi h L0=ln(a−e−b1⋅ 0)
Y=Y0⋅(1−e−b1⋅
1−e−b1⋅ 0)b2+b3/X0 (E4)
Ko
Y=a1⋅e−a2⋅ −a3
a
1
=X X0=Y0
e−b2⋅ −b3
0 Y=Y0⋅e−b2⋅ −b3
e−b2⋅ −b3
0
(E5)
a
2
=X X0= − ln(Y0/b1)/ −b3
0
Y=a1⋅(Y0
b1)( 0
)b3 (E6)
Weibull
Y=
a1⋅(1−e(− a2⋅ a3))
a
2
=X X0= − ln(1− (Y0/b1) )/ b3
1 Y=b1⋅(1− (1−Y0/b1))( / 0)b3 (E7)
Hoss eld
Y=a1
1+a2⋅ −a3
A
2
=X X0= −b3
0⋅(b1/Y0−1)Y=b1
1−(1−b1
Y0)⋅( 0
)b3
(E8)
Fig. 4. Resin yields o indi idual g oo es du ing he apping season o he wo
me hods, non-mechanized apping me hod (con inuous line) and mechanized
ci cula g oo e me hod (dashed line).
´
O. L´
opez-´
Al a ez e al.

Fo es Ecology and Managemen 549 (2023) 121501
6
om he minimum o he maximum p oduc ion in he pe iod . The
ep esen ed cu es a e hose o 187.00, 410.25, 633.50, 856.75, and
1080.00 accumula ed g ams a ime . Fig. 8 plo s he RRMSE ob ained
by he model when making p edic ions on es da a as a unc ion o ime.
The model achie es an RRMSE o less han 30 % o o ecas s up o 90
days (Fig. 8), which is e y good o o ecas s o less han one mon h,
good o o ecas s o wo mon hs and ai o o ecas s o h ee mon hs
acco ding o Li e al. (2013). This s a is ic ises abo e 30 % o o ecas s
longe han 90 days.
3.4. Mechanized ci cula g oo e me hod
The i s o he da a ob ained wi h he mechanized ci cula g oo e
me hod (Table 4), show a ou able esul s bu mo e modes compa ed o
he non-mechanized apping me hod. The R
2
alues achie ed wi h he
aining da a ange om 0.695 o 0.82. The RMSE alues and AIC alues
ange om 212.512 o 276.569 and om 16844.4 o 17496.84,
espec i ely. When assessing he models on es da a, he RRMSE was
again ob ained alues below 30 %. The bias a ied om 0.98 o 41.45 g,
and he MAE anged om 150.95 and 199.56 g.
Consis en ly, equa ions de i ed om he Be alan y-Richa ds model
demons a e supe io goodness-o - i s a is ics, yielding he mos
op imal i o he da a. Con e sely, he poo es i s a e obse ed wi h he
same equa ions as in he non-mechanized apping me hod case, speci -
ically he anamo phic and mul iple asymp o ic o m o he Be alan y-
Richa ds equa ion, as well as he Ko model by sol ing he “a
1
″
pa ame e .
The E4 o GADA o m o he Be alan y-Richa ds model is ound o
Fig. 5. Accumula ed esin yield du ing esin apping season o bo h esin apping me hods, non-mechanized apping me hod (g aphic a he op) and mechanized
ci cula g oo e me hod (g aphic a he bo om).
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
7
be he equa ion ha achie es he highes le el o s a is ical goodness o
i wi h he aining da a (Table 4). G aphical analysis o he esiduals
o he E4 model indica es he absence o any obse able he e o-
scedas ici y o au oco ela ion among he esiduals (Fig. 7). When es ed
agains he alida ion da a, he equa ion displayed he lowes RRMSE
alue o 29 %, a bias o 0.98 and a MAE o 150.95 g.
To ep esen he cu es o he adjus men s made wi h he equa ion
on he o iginal da a, i has been done in he same way as in he p e ious
Table 3
Pa ame e es ima es and goodness-o - i s a is ics o non-mechanized apping me hod.
Model T ain Tes
Pa ame e Es ima e S anda d e o R
2
RMSE AIC RRMSE (%) Bias MAE
E1 b
2
−0.003 0.001 0.818 295.901 18432.44 26.370 16.09 177.78
b
3
0.984 0.037
E2 b
1
4180.555 16.090 0.856 263.305 18129.92 24.806 −14.66 167.41
b
3
1.535 0.028
E3 b
1
15767.82 2.216x103 0.873 246.980 17964.02 23.189 −19.48 170.38
b
2
1.497x10−3 2.148x10−4
E4 b
1
1.608x10−3 1.425x10−4 0.873 247.330 17969.69 23.415 −24.34 170.85
b
2
−1.029x104 1.284x103
b
3
9.908x104 1.228x104
E5 b
2
17 3.582 0.812 300.758 18474.63 27.018 19.74 175.57
b
3
0.096 0.034
E6 b
1
3.447x104 7.237x103 0.857 261.883 18115.89 24.198 −11.32 164.38
b
3
0.318 0.019
E7 b
1
3.768x103 1.286x102 0.854 264.460 18141.26 25.047 −14.12 167.66
b
3
1.398 1.631x10−2
E8 b
1
5.233x103 2.449x102 0.855 263.686 18133.66 24.827 −14.19 167.72
b
3
1.465 2.202x10−2
Table 4
Pa ame e es ima es and goodness-o - i s a is ics o mechanized ci cula g oo e me hod.
Model T ain Tes
Pa ame e Es ima e S anda d e o R
2
RMSE AIC RRMSE (%) Bias MAE
E1 b
2
0.005 0.001 0.696 276.414 17495.44 41.971 41.45 199.43
b
3
1.247 0.059
E2 b
1
2.845x103 64.180 0.770 236.704 17110.2 33.698 12.67 170.69
b
3
1.589 0.029
E3 b
1
4.935x103 3.447x102 0.811 217.614 16901.33 29.844 4.29 156.77
b
2
3.603x10−3 3.076x10−4
E4 b
1
1.077x10−2 9.118x10−4 0.820 212.512 16844.4 29.244 0.98 150.95
b
2
−13.780 1.543
b
3
1.197x102 12.040
E5 b
2
11.931 0.629 0.695 276.569 17496.84 41.916 44.63 199.56
b
3
0.189 0.046
E6 b
1
5.885x103 4.574x102 0.817 214.433 16864.75 29.437 0.83 152.23
b
3
0.541 2.092x10−2
E7 b
1
2820.459 67.679 0.763 244.104 17186.67 35.540 21.69 178.76
b
3
1.357 0.017
E8 b
1
3.120x103 92.570 0.781 234.208 17083.87 33.020 12.54 167.14
b
3
1.563 0.024
Fig. 6. The le igu e was he age-dependen dynamic model E3 o non-mechanized no ching and he accumula ed esin yield cu es o 187.00, 410.25, 633.50,
856.75, and 1080.00 g a he ime o 56 days. The igh igu e is he Residuals s Residuals i +1 o non-mechanized no ching model E3.
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
8
case. In his case, he cu es ep esen ed a ime a e hose o 55.00,
333.75, 612.50, 891.25, and 1170.00 g (Fig. 7). The model ob ained an
RRMSE lowe han 35 % when he p edic ion lag was 85 days (Fig. 8),
also educing he mean p edic ion e o a he end o he apping season
below ha ob ained by he model i ed o he o he ex ac ion me hod.
In his case, o ecas s o less han one mon h could be desc ibed as good
and hose o less han wo and a hal mon hs as ai (Li e al., 2013;
Despo o ic e al., 2016).
4. Discussion
This s udy ep esen s he i s case o modelling indi idual ee
accumula i e esin yield h oughou he esin p oduc ion season. I
demons a es ha me hods simila o hose used o model di e en
na u al g ow h p ocesses can be success ully applied. The unp ece-
den ed na u e o his p oduc ion modelling poses a challenge o he
discussion, gi en he lack o p e ious s udies on esin o es ablish a
compa a i e amewo k. Fi ing models desc ibing he end o cumu-
la i e p oduc ion du ing he NTFP ha es season is challenging due o
he high a iabili y in indi idual ee p oduc ion (G´
omez-Ga cía e al.,
2022; L´
opez-Al a ez e al., 2023b) and he in luence ha en i onmen
changes ha e on he yield o his enewable p oduc (V´
azquez-Gonz´
alez
e al., 2021b).
Rega dless o he ex ac ion me hod, he equa ions de i ed om he
Be alan y-Richa ds model wi h ADA and GADA me hodologies we e
he ones ha p esen ed he bes pe o mance. Likewise, when
compa ing he s a is ics o he models based on he ex ac ion me hods,
he model o he non-mechanized me hod ob ained be e goodness-o -
i han he mechanized me hod.
The ob ained esul s sugges ha non-mechanized apping may lead
o a mo e p onounced yield inc ease in indi idual g oo es owa ds he
middle and end o he season when compa ed o he mechanized ci cula
g oo e me hod (Fig. 4). This beha io , de i ed mainly om he apping
me hod (L´
opez-´
Al a ez e al., 2023b), accen ua es he p esence o an
in lec ion poin in he cumula i e p oduc ion cu e, making i mo e
simila o g ow h p ocesses. In he non-mechanized apping me hod, he
equa ions ha ob ained he highes R
2
in he aining phase we e hose
de i ed om he Be alan y-Richa ds equa ion by ADA and GADA
me hodology (E3 and E4, espec i ely). The simila i y in desc ibing he
end o he da a be ween he ADA and GADA o ms o he Be alan y-
Richa ds model is common in o he o es species g ow h modelling
wo ks (T im e al., 2020; Wang e al., 2020). The GADA o m gene ally
ob ains be e goodness-o - i s a is ics, as i can gene a e polymo phic
cu es wi h mul iple asymp o es. This s a emen is e i ied in he case o
he non-mechanized ex ac ion me hod, as he GADA o mula ion ob-
ained he bes s a is ics, simila o he epo ed by S´
anchez-Gonz´
alez
Fig. 7. The le igu e was he age-dependen dynamic model E4 o mechanized ci cula g oo e and he accumula ed esin yield cu es o 55.00, 333.75, 612.50,
891.25, and 1170.00 g a he ime o 56 days. The igh igu e is he Residuals s Residuals i +1 o mechanized ci cula g oo e model E4.
Fig. 8. Rela i e oo mean squa e e o (RRMSE) by lag o p edic ion o mechanized ci cula g oo e (ci cles) and non-mechanized no ching ( iangles).
´
O. L´
opez-´
Al a ez e al.
Fo es Ecology and Managemen 549 (2023) 121501
9
e al. (2008) when modeling co k g ow h. The goodness-o - i s a is ics
ob ained in bo h i s a e simila o hose ob ained when using he GADA
o m o he Be alan y-Richa ds equa ion o model p oduc ion o e
ime o long o a ion species such as he peduncula e oak (G´
omez-Ga cía
e al., 2015), and sligh ly lowe han hose ob ained when modeling
sho o a ion species (Di´
eguez-A anda e al., 2005; Seki and Sakici,
2017) o co k p oduc ion (S´
anchez-Gonz´
alez e al., 2008).
The good pe o mance o he models in e ms o goodness-o - i , as
ob ained in his g ound-b eaking s udy, e i ies ha hese echniques
can be used o model esin p oduc ion du ing he season. To make hese
models a alid ool o modelling cumula i e esin p oduc ion, bo h in
he ield o ac i i y o o es wo ke s, in he esin apping indus y and in
he esin p ocessing indus y, i would be necessa y o de elop his
p elimina y s udy u he and o ca y ou a la ge sample, ep esen a-
i e o he whole en i onmen al g adien and including all he possible
casuis y unde which esin pine o es s can be ound. The e o e, based
on he abo e, desc ibing he indi idual ee accumula ed esin p oduc-
ion o P. pinas e can be included among he uses o he on Be alan y
model (Be alan y, 1949, 1957).
The p edic ions made wi h he es da a yielded a e age RRMSE
alues below 30 %, ega dless o he ex ac ion me hod. These a e age
alues can be classi ied as ai acco ding o Despo o ic e al. (2016) and
Li e al. (2013). Simpli ying and no conside ing he pa icula lag o
each p edic ion, he models ha e an a e age accu acy o 76.81 % and
70.76 %, o he non-mechanized and mechanized me hods, espec-
i ely. These accu acy pe cen ages a e sligh ly lowe han hose ach-
ie ed by he g ow h models used in dend ome y (P ada e al., 2019)
and a e mo e simila o hose ob ained by S´
anchez-Gonz´
alez e al.
(2008) modeling co k g ow h. When e alua ing he accu acy pe cen age
as a unc ion o he lag o he p edic ions, i was ound ha a p edic ion
made wi h an 80-day lag, which is usually hal o he esin- apping
campaign, had an e o a e o app oxima ely 35 %. These accu acy
pe cen ages as a unc ion o he lag a e sligh ly lowe han hose ach-
ie ed by dend ome ic models (Gea-Izquie do e al., 2008). The lowe
p ecision alues may be due o he small sample size o he high a i-
abili y o esin p oduc ion o e ime due o dasome ic and en i on-
men al a iables, as da a om dend ome ic a iables end o show less
a iabili y in p oduc ion be ween measu emen pe iods.
To acili a e scien i ic esea ch, one po en ial unc ionali y o he
esul s ob ained is hei use ulness o eco e missing da a om he
empo al weigh se ies. Collec ing in e media e weigh da a is a di icul
ask (Calama e al., 2011), and usually he e a e da a ha we e no
measu ed a he p ecise momen when each new s ip was done bu we e
accumula ed and weighed a e wo o mo e apping pe iods. In his
way, he accumula ed p oduc ion is accoun ed o , bu i is no known
which p oduc ion belongs o each pe iod. By using his ype o models,
he p oduc ion cu es can be ec ea ed, and he missing in e media e
da a can be eco e ed. They a e pa icula ly ele an o esea che s
who wish o gene a e and moni o p oduc ion cu es and ha e missing
da a in he p oduc ion se ies. The good p edic i e capaci y o he models
makes hem a alid ool o he mul i unc ional managemen o pine
o es s, which is easible as long as he owne s a e p o ided wi h he
app op ia e ools o ca y i ou (Sabas ian e al., 2019). By using hese
models, o es manage s will be able o p edic inpu s be o e he esin
ha es ing season ends, allowing hem o op imize he du a ion o he
campaign o maximise he annual income om esin sales. This makes
esin ex ac ion mo e a ac i e o o es p oduce s, c ea ing g een jobs
and p omo ing he bioeconomy by encou aging he p oduc ion o a
biop oduc such as pine esin.
5. Conclusion
Wi h he inc easing impo ance o esin p oduc ion in he NTFPs and
o es y sec o , he e is a need o u he de elop and expand ou
knowledge. This wo k demons a es o he i s ime ha age-in a ian
models a e a alid ool o he managemen o pine o es s unde esin
apping p ocedu es.
Acco ding o R
2
, RMSE, AIC, RRMSE, Bias and MAE, he ADA and
GADA o mula ions o he Be alan y-Richa ds model pe o med bes in
modelling cumula i e esin p oduc ion. This implies he inclusion o
accumula i e esin yield o indi idual ees as one o he many phe-
nomena ha can be modelled wi h hese equa ions. Fu he mo e, he
abili y o accu a ely p edic inal p oduc ion once he campaign has
s a ed using hese models ep esen s a signi ican ad ance in mul i-
unc ional o es managemen . In addi ion, his ype o model applied o
cumula i e annual esin p oduc ion can p o ide a basis o econ-
s uc ing missing da a se ies wi h a high deg ee o accu acy. To u he
ad ance his sec o , public o p i a e ini ia i es should p o ide suppo
and ollow he example o o es g ow h models by es ablishing a
ne wo k o pe manen plo s o pa ame ising na ional-scale esin p o-
duc ion a iables and de eloping mo e obus models. Enabling
ad ancemen in he e ec i e join managemen o wood esou ces and
esin p oduc ion, acili a ing he c ea ion o echnical–economic models
and he implemen a ion o inno a i e, cu ing-edge me hodologies in
cons uc ing said models.
Founding and acknowledgemen s
This wo k was suppo ed by he Spanish Go e nmen (“ACREMA”,
MAPA/AEI-Ag i/FEADER, UE) [O00000226e2000043659], he Gali-
cian Go e nmen (Xun a de Galicia) wi h a g an o Compe i i e
Re e ence G oups ED431C-2021-27 and he p e-doc o al con ac
Campus Te a-USC 2023. The au ho s hank FORESIN and CIF Lou izan,
o he ield assessmen s.
CRediT au ho ship con ibu ion s a emen
´
Osca L´
opez-´
Al a ez: Concep ualiza ion, Me hodology, Fo mal
analysis, W i ing – o iginal d a , Visualiza ion. Luis F anco-V´
azquez:
Me hodology, W i ing – e iew & edi ing. Manuel Ma ey-Pe ez:
Concep ualiza ion, Me hodology, Fo mal analysis, Resou ces, W i ing –
e iew & edi ing, Supe ision, P ojec adminis a ion, Funding
acquisi ion.
Decla a ion o Compe ing In e es
The au ho s decla e ha hey ha e no known compe ing inancial
in e es s o pe sonal ela ionships ha could ha e appea ed o in luence
he wo k epo ed in his pape .
Da a a ailabili y
Pe mission o sha e he da a mus be ob ained om he ins i u ions
ha unded he esea ch
Re e ences
Bailey, R.L., Clu e , J.L., 1974. Base-age in a ian polymo phic si e cu es. Fo . Sci. 20,
155–159. h ps://doi.o g/10.1093/ o es science/20.2.155.
Be alan y, L.V., 1957. Quan i a i e laws in me abolism and g ow h. Q. Re . Biol. 32,
217–231. h ps://doi.o g/10.1086/401873.
Cagla , M.U., Teu el, A.I., Wilke, C.O., 2018. Sicega : R package o sigmoidal and
double-sigmoidal cu e i ing. Pee J 6, e4251.
Calama R., Miina, J., de-Miguel, S., Bone , J. A., Mouni , F., Tom´
e, M., Ma ínez-
Jaú egui, M., He uzo, C., Pel ola, R., Salo, K., Ku ila M., He n´
andez-Rod íguez,
M., Ma ín-Pin o, P., S´
anchez-Gonz´
alez, M., 2020. Da a & models: impo ance o
assessing and o ecas ing non-wood o es p oduc s in Eu ope, in: Vacik, H., Hale,
M., Spiecke , H., Pe enella, D. Tom´
e, M., Non-Wood Fo es P oduc s in Eu ope.
Ecology and managemen o mush ooms, ee p oduc s, unde s o y plan s and
animal p oduc s. Ou comes o he COST Ac ion FP1203 on Eu opean NWFPs, BoD,
No de s ed , pp. 43-78.
Calama Sainz, R., Tome, M., S´
anchez-Gonz´
alez, M., Miina., J., Spanos, K., Palahi, M.,
2011. Modelling Non-Wood Fo es P oduc s in Eu ope: a e iew. Fo es Sys . 3 (4),
69.
´
O. L´
opez-´
Al a ez e al.