Learning to use ‘dyntaper’

Author

Affiliation

Roberto Scotti (scotti@uniss.it)

NuoroForestrySchool-DipAGR-UniSS.it

Published

12 / Aug / 2023

Introduzione

Oscar Garcia¹ ha teorizzato che la forma del fusto degli alberi sia modellizzabile in riferimento a proprietà ecofisiologiche da lui identificate².

Implementando tale teorizzazione l’autore ha pubblicato ‘dyntaper’, una libreria di funzioni R per produrre ed usare funzioni di profilo di fusto arborei.

Dettagli sulla teoria e sull’utilizzo del pacchetto sono presentati dall’autore in apposita vignette.

suppressPackageStartupMessages(library(tidyverse))
library(dyntaper)

Dati di esempio

Il pacchetto include ’brink' un piccolo data set con il rilievo del profilo di dieci fusti di diverse dimensioni.

# > str(brink)
# data.frame':  108 obs. of  5 variables:
#  $ Tree: Factor w/ 10 levels | Identificaivo del fusto
#  $ h   : num                 | altezza della sezione
#  $ dib : num                 | diametro (sotto corteccia) della sezione
#  $ Dob : num                 | diametro del fusto a petto d'uomo (sopra corteccia)
#  $ H   : num                 | altezza totale del fusto
#                              diametri in [cm], altezze in [m]
brink |> 
  summarise(num_sez = n(), .by = c(Tree, Dob, H))

   Tree  Dob     H num_sez
1   116 18.0 21.64      12
2   121 20.6 21.34      12
3   127 15.5 19.81      11
4   128 13.2 18.59      10
5   138 13.7 17.98      10
6   139 13.2 18.28      10
7   142 17.5 19.81      11
8   147 11.9 13.10       9
9   149 18.0 21.03      11
10  152 18.0 20.73      12

brink |> 
  ggplot(aes(h, dib, colour = Tree)) +
  geom_line() +geom_point()

Corteccia

La stima del profilo serve per valutare la possibibilità di utilizzare il toppo (la porzione di fusto) per produrre, in segheria, un determinato assortimento. Il diametro che interessa è quindi la misura ‘sotto-corteccia’.

Nel dataset di esempio le misure del diametro della sezione alle diverse altezze sono state rilevate sotto-corteccia. Come di consueto, la misura del diametro a petto d’uomo è rilevata sopra-corteccia.

Per focalizzare la stima dei profili sul diamtro sotto-corteccia, l’autore suggerisce di calcolare ‘k’ un fattore di riduzione medio da applicare poi a tutti i fusti considerati.

# calcolo fattore di riduzione del diametro a petto d'uomo 
# da sopra- a sotto- corteccia, utilizzando il rilievo sotto a 1.35 m
k <- brink |>
  filter(h == 1.35) |> 
  summarise(k = mean(dib/Dob)) |> 
  pluck(1)

cat(" *** k = ", k, "\n")

 *** k =  0.9081336

# in realtà l'autore, alla fine, sceglie di utilizzare un valore approssimato:
k <- 0.908

Stima dei coefficienti

La stima del valore dei parametri non è di semplice spiegazione, nel dettaglio.

Il procedimento suggerito nel documento di presentazione del pacchetto rappresenta un esempio molto semplificato, non direttamente replicabile in un contesto operativo. Tuttavia, essendo questa esposizione focalizzata sull’utilizzo delle funzioni e non sulla stima del valore dei coefficienti, l’esempio basta ed avanza.

#  *** Prima iterazione ***
(expexp <- nls(dib ~ taper(h, H, k*Dob, 
                          c(b1, 0, b3, b4, 0), 
                          1.35),
              data=brink, 
              start=c(b1=4, b3=1, b4=1)))

Nonlinear regression model
  model: dib ~ taper(h, H, k * Dob, c(b1, 0, b3, b4, 0), 1.35)
   data: brink
    b1     b3     b4 
2.0572 0.4311 0.6223 
 residual sum-of-squares: 33.77

Number of iterations to convergence: 10 
Achieved convergence tolerance: 2.908e-06

#  *** Interazione finale ***\n
#  utilizza e raffina le stime di prima iterazione, 
#  modificando il valore dei coefficienti inputati 'b2' e 'b5'
(hyplin <- nls(dib ~ taper(h, H, k*Dob, 
                           c(b1, -1, b3, b4, 1),
                           1.35), 
               data=brink, 
               start=coef(expexp)))

Nonlinear regression model
  model: dib ~ taper(h, H, k * Dob, c(b1, -1, b3, b4, 1), 1.35)
   data: brink
    b1     b3     b4 
2.2255 0.2826 1.0873 
 residual sum-of-squares: 31.42

Number of iterations to convergence: 5 
Achieved convergence tolerance: 2.953e-07

Composizione del vettore dei coefficienti ‘`b`’

coef(hyplin)

       b1        b3        b4 
2.2254819 0.2826045 1.0872791

(coeff <- c(coef(hyplin)["b1"], -1, coef(hyplin)["b3"], coef(hyplin)["b4"], 1))

        b1                    b3         b4            
 2.2254819 -1.0000000  0.2826045  1.0872791  1.0000000

Valutazione dell’accuratezza della stima dei diametri

brink |>
  mutate(dib_stimato = taper(h, H, Dob * k, b = coeff, bh = 1.3)) |>
  ggplot(aes(y = h)) +
  facet_grid(cols = vars(Tree), scales = "free_x") +
  xlab("d_sez [cm]") + ylab("h_sez [m]") +
  
  # misure
  geom_point(aes(x = dib),
             shape = 4,
             size = 1,
             colour = "blue") +
  geom_line(aes(x = dib), 
            linetype = "dotted", colour = "blue") +
  
  # diametri stimati
  geom_point(aes(x = dib_stimato),
             shape = 1,
             size = 2,
             colour = "grey") +
  geom_line(aes(x = dib_stimato),
            linetype = "solid",
            colour = "grey",
            linewidth = 1,
            alpha = .5) +
  ggtitle("Profili misurati (blue) e stimati (grey)")