# dendRoAnalyst: end-to-end workflow and function tour

#### Package vignette

## Overview

This vignette presents a complete workflow for dendRoAnalyst and also acts as a compact function tour for the user-facing functions covered by the uploaded source files. The sequence follows a typical dendrometer analysis pipeline:

1. read dendrometer and climate data,
2. check time resolution, jumps, and gaps,
3. resample and truncate,
4. smooth the signal,
5. compute daily summaries,
6. classify stem-cycle phases,
7. classify zero-growth phases,
8. perform superposed epoch analysis,
9. fit growth curves and compare methods,
10. detrend and standardize,
11. compute moving climate correlations,
12. analyze harsh-climate episodes,
13. interpolate missing values from a network, and
14. inspect periodicity with wavelet tools.

The vignette is written in a CRAN-friendly style. Lightweight sections are executed; heavier, interactive, or computationally expensive sections are shown as code examples with `eval = FALSE`. All the summary and plot functions can also be executed by simply applying `summary()` and `plot()` functions.

## Example datasets

    data(gf_nepa17)
    data(nepa17)
    data(ktm_clim_hourly)
    data(ktm_rain17)

    head(gf_nepa17[, 1:3])
    #>                  Time       T2       T3
    #> 1 2017-01-01 00:00:00 62.26459 48.16801
    #> 2 2017-01-01 01:00:00 62.27330 48.16559
    #> 3 2017-01-01 02:00:00 62.30913 48.16269
    #> 4 2017-01-01 03:00:00 62.34835 48.16414
    #> 5 2017-01-01 04:00:00 62.37740 48.16269
    #> 6 2017-01-01 05:00:00 62.42001 48.16705
    head(ktm_clim_hourly)
    #> # A tibble: 6 × 5
    #>   TIME                 temp  prec   VPD    RH
    #>   <dttm>              <dbl> <dbl> <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  7.02  6.28 0.162  83.8
    #> 2 2017-01-01 01:00:00  6.86  1.22 0.179  82.0
    #> 3 2017-01-01 02:00:00  6.80  2.41 0.154  84.4
    #> 4 2017-01-01 03:00:00  8.89  2.99 0.194  83.0
    #> 5 2017-01-01 04:00:00  9.96  3.61 0.282  77.0
    #> 6 2017-01-01 05:00:00 10.8   4.44 0.296  77.2

Throughout this vignette we use:

- `gf_nepa17` for subdaily dendrometer workflows,
- `nepa17[1:1000, ]` for compact daily examples,
- `ktm_clim_hourly` for climate-feature examples, and
- `ktm_rain17` for harsh-climate examples.

    dm_hourly <- gf_nepa17
    dm_daily_demo <- nepa17[1:1000, ]
    dm_phase_demo <- gf_nepa17[1:800, ]
    clim_hourly <- ktm_clim_hourly
    clim_rain <- ktm_rain17

## 1. Reading dendrometer and climate data

### `read.dendrometer()`

When users start from raw files, `read.dendrometer()` is the recommended entry point.

    raw_dm <- read.dendrometer(
      file = "inst/extdata/my_dendrometer_file.csv",
      detect_resolution = TRUE,
      return_report = TRUE,
      quiet = FALSE
    )

This call returns an imported dendrometer table and, optionally, an import report describing parsing and time-resolution diagnostics.

### `read.climate()`

    clim_std <- read.climate(
      x = clim_hourly,
      time_col = "TIME",
      vars = c("temp", "prec", "VPD", "RH"),
      verbose = FALSE
    )
    head(clim_std)
    #> # A tibble: 6 × 5
    #>   TIME                 temp  prec   VPD    RH
    #>   <dttm>              <dbl> <dbl> <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  7.02  6.28 0.162  83.8
    #> 2 2017-01-01 01:00:00  6.86  1.22 0.179  82.0
    #> 3 2017-01-01 02:00:00  6.80  2.41 0.154  84.4
    #> 4 2017-01-01 03:00:00  8.89  2.99 0.194  83.0
    #> 5 2017-01-01 04:00:00  9.96  3.61 0.282  77.0
    #> 6 2017-01-01 05:00:00 10.8   4.44 0.296  77.2

The standardized climate object can later be aggregated to daily scale, converted to subdaily climate features, or joined directly to dendrometer outputs.

## 2. Looking for jumps or gaps in dendrometer data

### `reso_dm()`

    reso_dm(dm_hourly[[1]])
    #> [1] 60

This returns the dominant sampling interval in minutes and helps diagnose irregular logging intervals before subsequent processing.

### `jump.locator()` and `i.jump.locator()`

    jump_free_manual <- jump.locator(
      df = dm_daily_demo[, 1:3],
      detection_method = "manual",
      manual_threshold = 1,
      adjustment_method = "window_median"
    )
    head(jump_free_manual)
    #>                  TIME       T2       T3
    #> 1 2017-01-01 00:00:00 62.26459 48.16801
    #> 2 2017-01-01 01:00:00 62.27330 48.16559
    #> 3 2017-01-01 02:00:00 62.30913 48.16269
    #> 4 2017-01-01 03:00:00 62.34835 48.16414
    #> 5 2017-01-01 04:00:00 62.37740 48.16269
    #> 6 2017-01-01 05:00:00 62.42001 48.16705

    jump_free_auto <- jump.locator(
      df = dm_daily_demo[, 1:3],
      detection_method = "auto",
      auto_method_penalty = 10,
      adjustment_method = "window_median"
    )

    # Interactive graphical workflow
    if (RUN_INTERACTIVE) {
      i.jump.locator(df = dm_daily_demo[, 1:3])
    }

### `dm.na.interpolation()` and its plot methods

    gap_demo <- dm_daily_demo[, 1:3]
    gap_demo[150:170, 2] <- NA

    gap_info <- dm.na.interpolation(
      df = gap_demo,
      resolution = 60,
      fill = FALSE,
      assess = FALSE
    )
    #> # A tibble: 1 × 9
    #>   type        series prev_time next_time start_missing
    #>   <chr>       <chr>  <dttm>    <dttm>    <dttm>
    #> 1 internal_na T2     NA        NA        2017-01-07 05:00:00
    #> # ℹ 4 more variables: end_missing <dttm>, n_consecutive <int>,
    #> #   duration_min <int>, duration_hr <dbl>
    head(gap_info$gap_info)
    #> # A tibble: 1 × 9
    #>   type        series prev_time next_time start_missing
    #>   <chr>       <chr>  <dttm>    <dttm>    <dttm>
    #> 1 internal_na T2     NA        NA        2017-01-07 05:00:00
    #> # ℹ 4 more variables: end_missing <dttm>, n_consecutive <int>,
    #> #   duration_min <int>, duration_hr <dbl>

    gap_fill <- dm.na.interpolation(
      df = gap_demo,
      resolution = 60,
      fill = TRUE,
      method = "spline",
      assess = TRUE,
      assess_lengths_hours = c(6, 12)
    )
    #> # A tibble: 1 × 9
    #>   type        series prev_time next_time start_missing
    #>   <chr>       <chr>  <dttm>    <dttm>    <dttm>
    #> 1 internal_na T2     NA        NA        2017-01-07 05:00:00
    #> # ℹ 4 more variables: end_missing <dttm>, n_consecutive <int>,
    #> #   duration_min <int>, duration_hr <dbl>
    head(gap_fill$data)
    #> # A tibble: 6 × 3
    #>   TIME                   T2    T3
    #>   <dttm>              <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  62.3  48.2
    #> 2 2017-01-01 01:00:00  62.3  48.2
    #> 3 2017-01-01 02:00:00  62.3  48.2
    #> 4 2017-01-01 03:00:00  62.3  48.2
    #> 5 2017-01-01 04:00:00  62.4  48.2
    #> 6 2017-01-01 05:00:00  62.4  48.2

    plot(gap_fill)

    plot(gap_fill)

    plot(gap_fill, original = gap_demo)

    plot(gap_fill)

These functions provide both the filled series and visual diagnostics of where and how the gaps were handled.

## 3. Resample and truncate the dataset

### `dendro.resample()`

    resampled_hourly <- dendro.resample(
      df = dm_hourly,
      by = "H",
      value = "mean",
      method = "auto"
    )
    head(resampled_hourly)
    #> # A tibble: 6 × 3
    #>   TIME                   T2    T3
    #>   <dttm>              <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  62.3  48.2
    #> 2 2017-01-01 01:00:00  62.3  48.2
    #> 3 2017-01-01 02:00:00  62.3  48.2
    #> 4 2017-01-01 03:00:00  62.3  48.2
    #> 5 2017-01-01 04:00:00  62.4  48.2
    #> 6 2017-01-01 05:00:00  62.4  48.2

### `dendro.truncate()`

    truncated_demo <- dendro.truncate(dm_hourly, CalYear = 2017, DOY = c(1, 120))
    head(truncated_demo)
    #> # A tibble: 6 × 3
    #>   TIME                   T2    T3
    #>   <dttm>              <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  62.3  48.2
    #> 2 2017-01-01 01:00:00  62.3  48.2
    #> 3 2017-01-01 02:00:00  62.3  48.2
    #> 4 2017-01-01 03:00:00  62.3  48.2
    #> 5 2017-01-01 04:00:00  62.4  48.2
    #> 6 2017-01-01 05:00:00  62.4  48.2

Resampling harmonizes time steps; truncation helps focus on one year or one season.

## 4. Smoothing

### `smooth_dm()`

    res_minutes <- reso_dm(dm_phase_demo[[1]])
    smoothed_medmean <- smooth_dm(
      time = dm_phase_demo[[1]],
      dm = dm_phase_demo[[2]],
      resolution_min = res_minutes,
      method = "median_mean",
      window_hours = 3
    )
    head(smoothed_medmean)
    #> [1] 62.26894 62.28234 62.31026 62.34496 62.38192 62.41565

    smoothed_pspline <- smooth_dm(
      time = dm_phase_demo[[1]],
      dm = dm_phase_demo[[2]],
      resolution_min = res_minutes,
      method = "pspline",
      window_hours = 6
    )

Smoothing is often applied before phase classification to suppress short spikes.

## 5. Daily approach and daily climate functions

### `daily.data()` and `plot.daily_output()`

    daily_out <- daily.data(df = dm_daily_demo, TreeNum = 1)
    head(daily_out)
    #> # A tibble: 6 × 14
    #>   DATE         Min Time_min   Max Time_max  mean median amplitude Time_min_h
    #>   <date>     <dbl> <chr>    <dbl> <chr>    <dbl>  <dbl>     <dbl>      <dbl>
    #> 1 2017-01-01  62.3 00:00:00  62.5 08:00:00  62.4   62.5    0.205           0
    #> 2 2017-01-02  62.5 14:00:00  62.5 08:00:00  62.5   62.5    0.0760         14
    #> 3 2017-01-03  62.5 00:00:00  62.8 15:00:00  62.6   62.6    0.243           0
    #> 4 2017-01-04  62.6 17:00:00  62.8 08:00:00  62.7   62.8    0.238          17
    #> 5 2017-01-05  62.6 15:00:00  62.8 09:00:00  62.7   62.6    0.149          15
    #> 6 2017-01-06  62.5 17:00:00  62.7 09:00:00  62.7   62.7    0.199          17
    #> # ℹ 5 more variables: Time_max_h <dbl>, lag_h <dbl>, Remarks <chr>,
    #> #   Max_diff <dbl>, Day_status <chr>

    plot(daily_out)

    plot(daily_out, type = "timing")

    plot(daily_out, type = "change")

    plot(daily_out, type = "timing_violin", Year = 2017, DOY = c(1, 6))
    plot(daily_out, type = "boxplot", stat = "amplitude", by = "month_of_year")
    plot(daily_out, type = "heatmap")

### `dm_daily_clim()` and `dm_join_daily_clim()`

    clim_day <- dm_daily_clim(
      clim_std,
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("VPD"),
      sum_vars  = c("prec"),
      lag_vars = c("VPD_max", "temp_mean"),
      lagmean_vars = c("temp_mean", "VPD_mean"),
      lagsum_vars = c("prec_sum"),
      lag_days = c(1, 3, 7)
    )
    head(clim_day)
    #> # A tibble: 6 × 21
    #>   DATE       temp_mean VPD_mean RH_mean VPD_max prec_sum VPD_max_lag_1d
    #>   <date>         <dbl>    <dbl>   <dbl>   <dbl>    <dbl>          <dbl>
    #> 1 2017-01-01      8.64    0.262    77.8   0.602  195.            NA
    #> 2 2017-01-02      8.70    0.358    69.9   0.711   45.2            0.602
    #> 3 2017-01-03      8.85    0.375    68.9   0.713    5.90           0.711
    #> 4 2017-01-04      9.55    0.474    63.4   0.879    0.236          0.713
    #> 5 2017-01-05      8.10    0.271    78.1   0.741   15.2            0.879
    #> 6 2017-01-06      8.12    0.190    84.3   0.478    3.39           0.741
    #> # ℹ 14 more variables: VPD_max_lag_3d <dbl>, VPD_max_lag_7d <dbl>,
    #> #   temp_mean_lag_1d <dbl>, temp_mean_lag_3d <dbl>, temp_mean_lag_7d <dbl>,
    #> #   temp_mean_lagmean_1d <dbl>, temp_mean_lagmean_3d <dbl>,
    #> #   temp_mean_lagmean_7d <dbl>, VPD_mean_lagmean_1d <dbl>,
    #> #   VPD_mean_lagmean_3d <dbl>, VPD_mean_lagmean_7d <dbl>,
    #> #   prec_sum_lagsum_1d <dbl>, prec_sum_lagsum_3d <dbl>,
    #> #   prec_sum_lagsum_7d <dbl>

    daily_clim_joined <- dm_join_daily_clim(daily_out, clim_day)
    head(daily_clim_joined)
    #> # A tibble: 6 × 34
    #>   DATE         Min Time_min   Max Time_max  mean median amplitude Time_min_h
    #>   <date>     <dbl> <chr>    <dbl> <chr>    <dbl>  <dbl>     <dbl>      <dbl>
    #> 1 2017-01-01  62.3 00:00:00  62.5 08:00:00  62.4   62.5    0.205           0
    #> 2 2017-01-02  62.5 14:00:00  62.5 08:00:00  62.5   62.5    0.0760         14
    #> 3 2017-01-03  62.5 00:00:00  62.8 15:00:00  62.6   62.6    0.243           0
    #> 4 2017-01-04  62.6 17:00:00  62.8 08:00:00  62.7   62.8    0.238          17
    #> 5 2017-01-05  62.6 15:00:00  62.8 09:00:00  62.7   62.6    0.149          15
    #> 6 2017-01-06  62.5 17:00:00  62.7 09:00:00  62.7   62.7    0.199          17
    #> # ℹ 25 more variables: Time_max_h <dbl>, lag_h <dbl>, Remarks <chr>,
    #> #   Max_diff <dbl>, Day_status <chr>, temp_mean <dbl>, VPD_mean <dbl>,
    #> #   RH_mean <dbl>, VPD_max <dbl>, prec_sum <dbl>, VPD_max_lag_1d <dbl>,
    #> #   VPD_max_lag_3d <dbl>, VPD_max_lag_7d <dbl>, temp_mean_lag_1d <dbl>,
    #> #   temp_mean_lag_3d <dbl>, temp_mean_lag_7d <dbl>, temp_mean_lagmean_1d <dbl>,
    #> #   temp_mean_lagmean_3d <dbl>, temp_mean_lagmean_7d <dbl>,
    #> #   VPD_mean_lagmean_1d <dbl>, VPD_mean_lagmean_3d <dbl>, …

### `dm_add_climate()` for daily outputs

    daily_clim_added <- dm_add_climate(
      daily_out,
      clim_hourly,
      scale = "daily",
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("temp", "VPD"),
      sum_vars  = c("prec")
    )
    head(daily_clim_added)
    #> # A tibble: 6 × 20
    #>   DATE         Min Time_min   Max Time_max  mean median amplitude Time_min_h
    #>   <date>     <dbl> <chr>    <dbl> <chr>    <dbl>  <dbl>     <dbl>      <dbl>
    #> 1 2017-01-01  62.3 00:00:00  62.5 08:00:00  62.4   62.5    0.205           0
    #> 2 2017-01-02  62.5 14:00:00  62.5 08:00:00  62.5   62.5    0.0760         14
    #> 3 2017-01-03  62.5 00:00:00  62.8 15:00:00  62.6   62.6    0.243           0
    #> 4 2017-01-04  62.6 17:00:00  62.8 08:00:00  62.7   62.8    0.238          17
    #> 5 2017-01-05  62.6 15:00:00  62.8 09:00:00  62.7   62.6    0.149          15
    #> 6 2017-01-06  62.5 17:00:00  62.7 09:00:00  62.7   62.7    0.199          17
    #> # ℹ 11 more variables: Time_max_h <dbl>, lag_h <dbl>, Remarks <chr>,
    #> #   Max_diff <dbl>, Day_status <chr>, temp_mean <dbl>, VPD_mean <dbl>,
    #> #   RH_mean <dbl>, temp_max <dbl>, VPD_max <dbl>, prec_sum <dbl>

## 6. Stem-cycle approach and associated functions

### `phase.sc()` and `plot.SC_output()`

    sc_out <- phase.sc(df = dm_phase_demo, TreeNum = 1)
    sc_out_sm <- phase.sc(df = dm_phase_demo, TreeNum = 1, smoothing = 12)
    head(sc_out$SC_phase)
    #> # A tibble: 6 × 3
    #>   TIME                   dm Phases
    #>   <dttm>              <dbl>  <dbl>
    #> 1 2017-01-01 00:00:00  62.3     NA
    #> 2 2017-01-01 01:00:00  62.3      3
    #> 3 2017-01-01 02:00:00  62.3      3
    #> 4 2017-01-01 03:00:00  62.3      3
    #> 5 2017-01-01 04:00:00  62.4      3
    #> 6 2017-01-01 05:00:00  62.4      3
    head(sc_out$SC_cycle)
    #> # A tibble: 6 × 8
    #>   Phases Start               End                 Duration_h Duration_m Magnitude
    #>    <dbl> <dttm>              <dttm>                   <dbl>      <dbl>     <dbl>
    #> 1      3 2017-01-01 00:00:00 2017-01-01 08:00:00          8        480   0.205
    #> 2      1 2017-01-01 08:00:00 2017-01-01 12:00:00          4        240  -0.0223
    #> 3      2 2017-01-01 12:00:00 2017-01-01 15:00:00          3        180   0.00871
    #> 4      1 2017-01-01 15:00:00 2017-01-01 16:00:00          1         60  -0.00290
    #> 5      2 2017-01-01 16:00:00 2017-01-01 20:00:00          4        240   0.0121
    #> 6      1 2017-01-01 20:00:00 2017-01-01 21:00:00          1         60  -0.00145
    #> # ℹ 2 more variables: rate <dbl>, DOY <dbl>

    plot(sc_out)

    plot(sc_out, type = "ribbon")

    plot(sc_out, type = "increment")

    plot(sc_out, type = "transition", x_axis = "doy", Year = 2017, DOY = c(50, 100))
    plot(sc_out, type = "balance", temporal = "daily")
    plot(sc_out, type = "frequency", temporal = "monthly")
    plot(sc_out, type = "heatmap", temporal = "daily")
    plot(sc_out, type = "boxplot", stat = "Magnitude", temporal = "monthly")

### `dm_subdaily_clim()`, `dm_join_subdaily_clim()`, `dm_join_phase_clim()`, and `dm_add_climate()`

    clim_sub <- dm_subdaily_clim(
      clim_hourly,
      mean_vars = c("temp", "VPD", "RH"),
      sum_vars  = c("prec"),
      lag_vars  = c("temp", "VPD", "RH"),
      roll_hours = c(1, 3, 6, 24),
      lag_hours  = c(1, 3, 6, 24)
    )
    head(clim_sub)
    #> # A tibble: 6 × 33
    #>   TIME                 temp  prec   VPD    RH temp_rollmean_1h temp_rollmean_3h
    #>   <dttm>              <dbl> <dbl> <dbl> <dbl>            <dbl>            <dbl>
    #> 1 2017-01-01 00:00:00  7.02  6.28 0.162  83.8             7.02            NA
    #> 2 2017-01-01 01:00:00  6.86  1.22 0.179  82.0             6.86            NA
    #> 3 2017-01-01 02:00:00  6.80  2.41 0.154  84.4             6.80             6.89
    #> 4 2017-01-01 03:00:00  8.89  2.99 0.194  83.0             8.89             7.52
    #> 5 2017-01-01 04:00:00  9.96  3.61 0.282  77.0             9.96             8.55
    #> 6 2017-01-01 05:00:00 10.8   4.44 0.296  77.2            10.8              9.88
    #> # ℹ 26 more variables: temp_rollmean_6h <dbl>, temp_rollmean_24h <dbl>,
    #> #   VPD_rollmean_1h <dbl>, VPD_rollmean_3h <dbl>, VPD_rollmean_6h <dbl>,
    #> #   VPD_rollmean_24h <dbl>, RH_rollmean_1h <dbl>, RH_rollmean_3h <dbl>,
    #> #   RH_rollmean_6h <dbl>, RH_rollmean_24h <dbl>, prec_rollsum_1h <dbl>,
    #> #   prec_rollsum_3h <dbl>, prec_rollsum_6h <dbl>, prec_rollsum_24h <dbl>,
    #> #   temp_lag_1h <dbl>, temp_lag_3h <dbl>, temp_lag_6h <dbl>,
    #> #   temp_lag_24h <dbl>, VPD_lag_1h <dbl>, VPD_lag_3h <dbl>, VPD_lag_6h <dbl>, …

    sc_phase_clim <- dm_join_phase_clim(
      sc_out,
      clim_hourly,
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("temp", "VPD"),
      sum_vars  = c("prec")
    )
    head(sc_phase_clim$SC_cycle)
    #> # A tibble: 6 × 15
    #>   Phases Start               End                 Duration_h Duration_m Magnitude
    #>    <dbl> <dttm>              <dttm>                   <dbl>      <dbl>     <dbl>
    #> 1      3 2017-01-01 00:00:00 2017-01-01 08:00:00          8        480   0.205
    #> 2      1 2017-01-01 08:00:00 2017-01-01 12:00:00          4        240  -0.0223
    #> 3      2 2017-01-01 12:00:00 2017-01-01 15:00:00          3        180   0.00871
    #> 4      1 2017-01-01 15:00:00 2017-01-01 16:00:00          1         60  -0.00290
    #> 5      2 2017-01-01 16:00:00 2017-01-01 20:00:00          4        240   0.0121
    #> 6      1 2017-01-01 20:00:00 2017-01-01 21:00:00          1         60  -0.00145
    #> # ℹ 9 more variables: rate <dbl>, DOY <dbl>, n_clim <int>,
    #> #   temp_mean_phase <dbl>, VPD_mean_phase <dbl>, RH_mean_phase <dbl>,
    #> #   temp_max_phase <dbl>, VPD_max_phase <dbl>, prec_sum_phase <dbl>

    sc_point_clim <- dm_join_subdaily_clim(sc_out_sm, clim_sub)
    head(sc_point_clim$SC_phase)
    #> # A tibble: 6 × 35
    #>   TIME                   dm Phases  temp  prec   VPD    RH temp_rollmean_1h
    #>   <dttm>              <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl>            <dbl>
    #> 1 2017-01-01 00:00:00  62.3     NA  7.02  6.28 0.162  83.8             7.02
    #> 2 2017-01-01 01:00:00  62.3      3  6.86  1.22 0.179  82.0             6.86
    #> 3 2017-01-01 02:00:00  62.3      3  6.80  2.41 0.154  84.4             6.80
    #> 4 2017-01-01 03:00:00  62.3      3  8.89  2.99 0.194  83.0             8.89
    #> 5 2017-01-01 04:00:00  62.4      3  9.96  3.61 0.282  77.0             9.96
    #> 6 2017-01-01 05:00:00  62.4      3 10.8   4.44 0.296  77.2            10.8
    #> # ℹ 27 more variables: temp_rollmean_3h <dbl>, temp_rollmean_6h <dbl>,
    #> #   temp_rollmean_24h <dbl>, VPD_rollmean_1h <dbl>, VPD_rollmean_3h <dbl>,
    #> #   VPD_rollmean_6h <dbl>, VPD_rollmean_24h <dbl>, RH_rollmean_1h <dbl>,
    #> #   RH_rollmean_3h <dbl>, RH_rollmean_6h <dbl>, RH_rollmean_24h <dbl>,
    #> #   prec_rollsum_1h <dbl>, prec_rollsum_3h <dbl>, prec_rollsum_6h <dbl>,
    #> #   prec_rollsum_24h <dbl>, temp_lag_1h <dbl>, temp_lag_3h <dbl>,
    #> #   temp_lag_6h <dbl>, temp_lag_24h <dbl>, VPD_lag_1h <dbl>, …

    sc_added_clim <- dm_add_climate(
      sc_out_sm,
      clim_hourly,
      scale = "subdaily",
      sub_mean_vars = c("temp", "VPD", "RH"),
      sub_sum_vars  = c("prec"),
      sub_lag_vars  = c("temp", "VPD", "RH"),
      roll_hours = c(1, 3, 6, 24),
      lag_hours  = c(1, 3, 6, 24)
    )
    head(sc_added_clim$SC_phase)
    #> # A tibble: 6 × 35
    #>   TIME                   dm Phases  temp  prec   VPD    RH temp_rollmean_1h
    #>   <dttm>              <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl>            <dbl>
    #> 1 2017-01-01 00:00:00  62.3     NA  7.02  6.28 0.162  83.8             7.02
    #> 2 2017-01-01 01:00:00  62.3      3  6.86  1.22 0.179  82.0             6.86
    #> 3 2017-01-01 02:00:00  62.3      3  6.80  2.41 0.154  84.4             6.80
    #> 4 2017-01-01 03:00:00  62.3      3  8.89  2.99 0.194  83.0             8.89
    #> 5 2017-01-01 04:00:00  62.4      3  9.96  3.61 0.282  77.0             9.96
    #> 6 2017-01-01 05:00:00  62.4      3 10.8   4.44 0.296  77.2            10.8
    #> # ℹ 27 more variables: temp_rollmean_3h <dbl>, temp_rollmean_6h <dbl>,
    #> #   temp_rollmean_24h <dbl>, VPD_rollmean_1h <dbl>, VPD_rollmean_3h <dbl>,
    #> #   VPD_rollmean_6h <dbl>, VPD_rollmean_24h <dbl>, RH_rollmean_1h <dbl>,
    #> #   RH_rollmean_3h <dbl>, RH_rollmean_6h <dbl>, RH_rollmean_24h <dbl>,
    #> #   prec_rollsum_1h <dbl>, prec_rollsum_3h <dbl>, prec_rollsum_6h <dbl>,
    #> #   prec_rollsum_24h <dbl>, temp_lag_1h <dbl>, temp_lag_3h <dbl>,
    #> #   temp_lag_6h <dbl>, temp_lag_24h <dbl>, VPD_lag_1h <dbl>, …

## 7. Zero-growth approach and associated functions

### `phase.zg()`, `twd.maxima()`, and `plot.ZG_output()`

    zg_out <- phase.zg(df = dm_phase_demo, TreeNum = 1)
    zg_out_sm <- phase.zg(df = dm_phase_demo, TreeNum = 1, smoothing = 6)
    head(zg_out$ZG_phase)
    #> # A tibble: 6 × 5
    #>   TIME                   dm Phases   TWD   GRO
    #>   <dttm>              <dbl>  <int> <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  62.3      2     0  62.3
    #> 2 2017-01-01 01:00:00  62.3      2     0  62.3
    #> 3 2017-01-01 02:00:00  62.3      2     0  62.3
    #> 4 2017-01-01 03:00:00  62.3      2     0  62.3
    #> 5 2017-01-01 04:00:00  62.4      2     0  62.4
    #> 6 2017-01-01 05:00:00  62.4      2     0  62.4
    head(zg_out$ZG_cycle)
    #> # A tibble: 6 × 14
    #>   Phases Start               End                 Duration_h Magnitude  rate
    #>    <int> <dttm>              <dttm>                   <dbl>     <dbl> <dbl>
    #> 1      2 2017-01-01 00:00:00 2017-01-01 09:00:00          9    0.205  22.8
    #> 2      1 2017-01-01 09:00:00 2017-01-02 00:00:00         15   NA      NA
    #> 3      2 2017-01-02 00:00:00 2017-01-02 09:00:00          9    0.0731  8.12
    #> 4      1 2017-01-02 09:00:00 2017-01-03 02:00:00         17   NA      NA
    #> 5      2 2017-01-03 02:00:00 2017-01-03 10:00:00          8    0.0692  8.65
    #> 6      1 2017-01-03 10:00:00 2017-01-03 14:00:00          4   NA      NA
    #> # ℹ 8 more variables: max.twd <dbl>, Max.twd.time <dttm>, AUC.load <dbl>,
    #> #   AUC.total <dbl>, ABr.value <dbl>, Avg.twd <dbl>, STD.twd <dbl>, DOY <dbl>

    max_twd_tbl <- twd.maxima(dm_phase_demo, TreeNum = 1, smoothing = 5)
    head(max_twd_tbl)
    #> # A tibble: 6 × 7
    #>   Start.time          End.time            TIME                Phases    TWD
    #>   <dttm>              <dttm>              <dttm>               <int>  <dbl>
    #> 1 2017-01-01 09:00:00 2017-01-02 00:00:00 2017-01-01 12:00:00      1 0.0223
    #> 2 2017-01-02 09:00:00 2017-01-03 02:00:00 2017-01-02 14:00:00      1 0.0760
    #> 3 2017-01-03 10:00:00 2017-01-03 14:00:00 2017-01-03 13:00:00      1 0.0692
    #> 4 2017-01-04 09:00:00 2017-01-19 07:00:00 2017-01-04 16:00:00      1 0.236
    #> 5 2017-01-04 09:00:00 2017-01-19 07:00:00 2017-01-04 17:00:00      1 0.238
    #> 6 2017-01-04 09:00:00 2017-01-19 07:00:00 2017-01-05 15:00:00      1 0.216
    #> # ℹ 2 more variables: Duration_from_start <drtn>, twd.number <int>

    plot(zg_out)

    plot(zg_out, temporal = "daily")

    plot(zg_out, temporal = "monthly", type = "twd")

    plot(zg_out, temporal = "monthly", type = "abr")
    plot(zg_out, type = "transition", x_axis = "doy", Year = 2017, DOY = c(50, 100))
    plot(zg_out, type = "boxplot", stat = "max.twd", temporal = "daily")

### Climate joins and climate-phase plots

    zg_phase_clim <- dm_join_phase_clim(
      zg_out,
      clim_hourly,
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("temp", "VPD"),
      sum_vars  = c("prec")
    )
    head(zg_phase_clim$ZG_cycle)
    #> # A tibble: 6 × 21
    #>   Phases Start               End                 Duration_h Magnitude  rate
    #>    <int> <dttm>              <dttm>                   <dbl>     <dbl> <dbl>
    #> 1      2 2017-01-01 00:00:00 2017-01-01 09:00:00          9    0.205  22.8
    #> 2      1 2017-01-01 09:00:00 2017-01-02 00:00:00         15   NA      NA
    #> 3      2 2017-01-02 00:00:00 2017-01-02 09:00:00          9    0.0731  8.12
    #> 4      1 2017-01-02 09:00:00 2017-01-03 02:00:00         17   NA      NA
    #> 5      2 2017-01-03 02:00:00 2017-01-03 10:00:00          8    0.0692  8.65
    #> 6      1 2017-01-03 10:00:00 2017-01-03 14:00:00          4   NA      NA
    #> # ℹ 15 more variables: max.twd <dbl>, Max.twd.time <dttm>, AUC.load <dbl>,
    #> #   AUC.total <dbl>, ABr.value <dbl>, Avg.twd <dbl>, STD.twd <dbl>, DOY <dbl>,
    #> #   n_clim <int>, temp_mean_phase <dbl>, VPD_mean_phase <dbl>,
    #> #   RH_mean_phase <dbl>, temp_max_phase <dbl>, VPD_max_phase <dbl>,
    #> #   prec_sum_phase <dbl>

    zg_point_clim <- dm_join_subdaily_clim(zg_out_sm, clim_sub)
    head(zg_point_clim$ZG_phase)
    #> # A tibble: 6 × 37
    #>   TIME                   dm Phases   TWD   GRO  temp  prec   VPD    RH
    #>   <dttm>              <dbl>  <int> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
    #> 1 2017-01-01 00:00:00  62.3      2     0  62.3  7.02  6.28 0.162  83.8
    #> 2 2017-01-01 01:00:00  62.3      2     0  62.3  6.86  1.22 0.179  82.0
    #> 3 2017-01-01 02:00:00  62.3      2     0  62.3  6.80  2.41 0.154  84.4
    #> 4 2017-01-01 03:00:00  62.3      2     0  62.3  8.89  2.99 0.194  83.0
    #> 5 2017-01-01 04:00:00  62.4      2     0  62.4  9.96  3.61 0.282  77.0
    #> 6 2017-01-01 05:00:00  62.4      2     0  62.4 10.8   4.44 0.296  77.2
    #> # ℹ 28 more variables: temp_rollmean_1h <dbl>, temp_rollmean_3h <dbl>,
    #> #   temp_rollmean_6h <dbl>, temp_rollmean_24h <dbl>, VPD_rollmean_1h <dbl>,
    #> #   VPD_rollmean_3h <dbl>, VPD_rollmean_6h <dbl>, VPD_rollmean_24h <dbl>,
    #> #   RH_rollmean_1h <dbl>, RH_rollmean_3h <dbl>, RH_rollmean_6h <dbl>,
    #> #   RH_rollmean_24h <dbl>, prec_rollsum_1h <dbl>, prec_rollsum_3h <dbl>,
    #> #   prec_rollsum_6h <dbl>, prec_rollsum_24h <dbl>, temp_lag_1h <dbl>,
    #> #   temp_lag_3h <dbl>, temp_lag_6h <dbl>, temp_lag_24h <dbl>, …

    zg_phase_added <- dm_add_climate(
      zg_out,
      clim_hourly,
      scale = "phase",
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("temp", "VPD"),
      sum_vars  = c("prec")
    )
    head(zg_phase_added$ZG_cycle)
    #> # A tibble: 6 × 21
    #>   Phases Start               End                 Duration_h Magnitude  rate
    #>    <int> <dttm>              <dttm>                   <dbl>     <dbl> <dbl>
    #> 1      2 2017-01-01 00:00:00 2017-01-01 09:00:00          9    0.205  22.8
    #> 2      1 2017-01-01 09:00:00 2017-01-02 00:00:00         15   NA      NA
    #> 3      2 2017-01-02 00:00:00 2017-01-02 09:00:00          9    0.0731  8.12
    #> 4      1 2017-01-02 09:00:00 2017-01-03 02:00:00         17   NA      NA
    #> 5      2 2017-01-03 02:00:00 2017-01-03 10:00:00          8    0.0692  8.65
    #> 6      1 2017-01-03 10:00:00 2017-01-03 14:00:00          4   NA      NA
    #> # ℹ 15 more variables: max.twd <dbl>, Max.twd.time <dttm>, AUC.load <dbl>,
    #> #   AUC.total <dbl>, ABr.value <dbl>, Avg.twd <dbl>, STD.twd <dbl>, DOY <dbl>,
    #> #   n_clim <int>, temp_mean_phase <dbl>, VPD_mean_phase <dbl>,
    #> #   RH_mean_phase <dbl>, temp_max_phase <dbl>, VPD_max_phase <dbl>,
    #> #   prec_sum_phase <dbl>

    plot(zg_phase_added, climate_var = "VPD_mean_phase", scale = "cycle", type = "boxplot", temporal = 'daily')

    plot(zg_phase_added, climate_var = "temp_mean_phase", scale = "cycle", type = "timeseries")

    plot(
      zg_phase_added,
      climate_vars = c("temp_mean_phase", "VPD_mean_phase", "RH_mean_phase"),
      scale = "cycle",
      type = "boxplot"
    )

    plot(
      zg_phase_added,
      scale = "cycle",
      type = "cor_heatmap",
      numeric_vars = c("Duration_h", "max.twd", "ABr.value",
                       "temp_mean_phase", "VPD_mean_phase", "RH_mean_phase")
    )

## 8. Superposed epoch and associated functions

### `dm_event_times()`, `dm_epoch_extract()`, `dm_epoch_test()`, and method dispatch

    event_tbl <- dm_event_times(zg_out, event = "all")
    head(event_tbl)
    #> # A tibble: 6 × 19
    #>   event_id event_time          event_type phase Phases Start
    #>      <int> <dttm>              <chr>      <chr>  <int> <dttm>
    #> 1        1 2017-01-01 00:00:00 GRO_start  GRO        2 2017-01-01 00:00:00
    #> 2        2 2017-01-01 09:00:00 TWD_start  TWD        1 2017-01-01 09:00:00
    #> 3        3 2017-01-01 09:00:00 GRO_end    GRO        2 2017-01-01 00:00:00
    #> 4        4 2017-01-01 12:00:00 MaxTWD     TWD        1 2017-01-01 09:00:00
    #> 5        5 2017-01-02 00:00:00 TWD_end    TWD        1 2017-01-01 09:00:00
    #> 6        6 2017-01-02 00:00:00 GRO_start  GRO        2 2017-01-02 00:00:00
    #> # ℹ 13 more variables: End <dttm>, Duration_h <dbl>, Magnitude <dbl>,
    #> #   rate <dbl>, max.twd <dbl>, Max.twd.time <dttm>, AUC.load <dbl>,
    #> #   AUC.total <dbl>, ABr.value <dbl>, Avg.twd <dbl>, STD.twd <dbl>, DOY <dbl>,
    #> #   phase_label <chr>

    epoch_extract <- dm_epoch_extract(
      x = zg_out,
      clim = clim_hourly,
      event = "phase_start",
      phase = "all",
      response_var = c("temp", "VPD"),
      before = 24,
      after = 24,
      unit = "hours"
    )
    head(epoch_extract)

    epoch_test <- dm_epoch_test(
      x = zg_out,
      clim = clim_hourly,
      event = "phase_start",
      phase = "all",
      response_var = c("temp", "VPD"),
      before = 24,
      after = 24,
      unit = "hours",
      n_boot = 99
    )

    epoch_summary <- summary(epoch_test)   # summary.dm_epoch()
    print(epoch_summary)                   # print.summary.dm_epoch()
    plot(epoch_test)

## 9. Event-based climate functions

    evt_zg <- dm_event_climate(
      zg_out,
      clim_hourly,
      event = "phase_start",
      phase = "all",
      windows = c(0, 3, 6, 12, 24),
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("VPD"),
      sum_vars  = c("prec")
    )
    head(evt_zg)

    evt_summary <- dm_event_climate_summary(
      evt_zg,
      climate_var = "VPD_mean_prev_6h",
      group_var = "Phase",
      by = "month_of_year"
    )
    head(evt_summary)

    evt_test <- dm_event_climate_test(
      evt_zg,
      climate_var = "VPD_mean_prev_6h",
      group_var = "Phase",
      by = "month_of_year",
      method = "auto"
    )
    head(evt_test$summary)
    head(evt_test$overall_tests)
    head(evt_test$pairwise_tests)

    plot(
      evt_zg,
      climate_var = "VPD_mean_prev_6h",
      group_var = "Phase",
      facet_by = "month_of_year",
      geom = "both",
      add_test = TRUE
    )

    evt_maxtwd <- dm_event_climate(
      zg_out,
      clim_hourly,
      event = "max_twd",
      windows = c(0, 3, 6, 12, 24),
      mean_vars = c("temp", "VPD", "RH"),
      max_vars  = c("VPD"),
      sum_vars  = c("prec")
    )

    plot(
      evt_maxtwd,
      climate_var = "VPD_mean_prev_6h",
      response_var = "max.twd",
      group_var = "Phase"
    )

## 10. Growth fit and evaluation

### `dm.fit.gompertz()`

    gomp_legacy <- dm.fit.gompertz(df = dm_hourly, TreeNum = 1, CalYear = 2017, verbose = FALSE)
    head(gomp_legacy)

### `dm.growth.fit()`, `plot.dm_growth_fit()`, and print/summary methods

    fit_single <- dm.growth.fit(
      df = dm_hourly,
      TreeNum = 1:2,
      method = "gompertz",
      year_mode = "yearly",
      verbose = FALSE
    )

    print(fit_single)                 # print.dm_growth_fit()
    fit_single_summary <- summary(fit_single)   # summary.dm_growth_fit()
    print(fit_single_summary)         # print.summary.dm_growth_fit()

    plot(fit_single, type = "fit")
    plot(fit_single, type = "season")
    plot(fit_single, type = "residuals")
    plot(fit_single, type = "timing")
    plot(fit_single, type = "overlay")
    plot(fit_single, type = "parameters")

### `dm.growth.fit.double()` and `dm.growth.evaluate()`

    fit_double <- dm.growth.fit.double(
      df = dm_hourly,
      TreeNum = 1:2,
      method = "gompertz",
      year_mode = "yearly",
      verbose = FALSE
    )
    print(fit_double)

    growth_eval <- dm.growth.evaluate(
      df = dm_hourly,
      TreeNum = 1:2,
      methods = c("gompertz", "loess", "spline", "double_gompertz"),
      year_mode = "yearly",
      verbose = FALSE
    )
    head(growth_eval)

    plot(growth_eval, metric = "rmse", type = "boxplot")
    plot(growth_eval, metric = "r2", type = "heatmap")

## 11. Detrending and associated functions

    detrended_fit <- dm.detrend.fit(fit_single)
    plot(detrended_fit, type = "compare")
    plot(detrended_fit, type = "fit")
    plot(detrended_fit, type = "residual")
    plot(detrended_fit, type = "detrended")
    plot(detrended_fit, type = "boxplot")

    mean_det <- mean_detrended.dm(detrended_fit)
    head(mean_det)
    plot(mean_det)

### `dm_standardize()` and `plot.dm_standardized()`

    std_nh <- dm_standardize(
      df = dm_hourly,
      season_type = "NH",
      method = "robust_amplitude"
    )
    head(std_nh$data)
    #> # A tibble: 6 × 5
    #>   TIME                season_year in_season       T2         T3
    #>   <dttm>                    <dbl> <lgl>        <dbl>      <dbl>
    #> 1 2017-01-01 00:00:00        2017 TRUE      0         0
    #> 2 2017-01-01 01:00:00        2017 TRUE      0.000566 -0.000200
    #> 3 2017-01-01 02:00:00        2017 TRUE      0.00289  -0.000439
    #> 4 2017-01-01 03:00:00        2017 TRUE      0.00544  -0.000319
    #> 5 2017-01-01 04:00:00        2017 TRUE      0.00733  -0.000439
    #> 6 2017-01-01 05:00:00        2017 TRUE      0.0101   -0.0000798
    head(std_nh$parameters)
    #> # A tibble: 2 × 17
    #>   tree  season_year season_type method        ref_type ref_n_days ref_doys_start
    #>   <chr>       <dbl> <chr>       <chr>         <chr>         <int>          <int>
    #> 1 T2           2017 NH          robust_ampli… first_v…         NA             NA
    #> 2 T3           2017 NH          robust_ampli… first_v…         NA             NA
    #> # ℹ 10 more variables: ref_doys_end <int>, q_low <dbl>, q_high <dbl>,
    #> #   n_total <int>, n_nonmiss <int>, season_start <dttm>, season_end <dttm>,
    #> #   reference_value <dbl>, center_value <dbl>, denominator <dbl>

    plot(std_nh)

## 12. Moving correlations and associated functions

    mov_out <- mov.cor.dm(
      df = dm_hourly,
      Clim = clim_hourly,
      TreeNum = 1,
      win_size = 15,
      cor_method = "pearson",
      boot = FALSE
    )

    print(mov_out)                # print.mov_cor_dm()
    mov_summary <- summary(mov_out)   # summary.mov_cor_dm()
    print(mov_summary)            # print.summary_mov_cor_dm()

    plot(mov_out, type = "heatmap")
    plot(mov_out, type = "line", annotate_peak = TRUE)
    plot(mov_out, type = "facet")
    plot(mov_out, type = "peak")
    plot(mov_out, type = "line")
    plot(mov_summary, type = "peak_corr")

    mov_boot <- mov.cor.dm(
      df = dm_hourly,
      Clim = clim_hourly,
      TreeNum = 1,
      win_size = 15,
      cor_method = "pearson",
      boot = TRUE,
      R = 99,
      set_seed = 1
    )
    plot(mov_boot, type = "line", show_ci = TRUE)

## 13. Harsh climate using `clim.twd()` and associated functions

    rel_out <- clim.twd(
      df = dm_hourly,
      Clim = clim_rain,
      thresholdClim = "<10",
      thresholdDays = ">5",
      showPlot = FALSE
    )

    clim_stats <- clim.twd.stats(
      rel_out,
      response = "full_period_change",
      group_by = "tree"
    )

    clim_stats_summary <- summary(clim_stats)   # summary.clim_twd_stats()
    print(clim_stats_summary)                   # print.summary.clim_twd_stats()
    plot(clim_stats, type = "trajectory")

    tree_info <- data.frame(
      tree = c("T2", "T3"),
      species = c("Pinus", "Pinus"),
      site = c("Ktm", "Ktm")
    )

    clim_test <- clim.twd.test(
      rel_out,
      tree_info = tree_info,
      parameter = "lag_to_return_zero",
      compare_by = "species",
      within = "year"
    )

    clim_test_summary <- summary(clim_test)   # summary.clim_twd_test()
    print(clim_test_summary)                  # print.summary.clim_twd_test()
    plot(clim_test)

## 14. Network interpolation

### `network.interpolation()` and `plot.network_interpolation()`

    df1 <- dm_hourly
    df1[40:50, "T2"] <- NA

    ref_net <- cbind(dm_hourly, dm_hourly[, 2:3], dm_hourly[, 2:3])
    colnames(ref_net) <- c("Time", "T1", "T2", "T3", "T4", "T5", "T6")

    net_interp <- network.interpolation(
      df1,
      ref_net,
      niMethod = "proportional",
      n_boot = 100,
      return_flags = TRUE,
      return_pi = TRUE,
      assess = TRUE,
      assess_lengths_steps = c(1, 2, 4),
      correct_gap_jumps = TRUE,
      jump_threshold = 0.05
    )
    head(net_interp)
    # to get the assessment stats we can use attr() function
    attr(out, "network_validation_summary")
    attr(out, "network_validation_points")
    attr(out, "network_validation_seasonal")

    plot(net_interp)
    plot(net_interp, type = "uncertainty")
    plot(net_interp, type = "availability")
    plot(net_interp, type = "summary", summary_metric = "n_gap_jumps_removed")
    plot(net_interp, type = "validation", validation_metric = "MAPE")
    plot(net_interp, type = "coverage")
    plot(net_interp, type = "compare")
    plot(net_interp, type = "seasonal_error", seasonal_metric = "mean_abs_error")

## 15. Wavelet analysis and associated functions

### `dm_wavelet()`, `summary.dm_wavelet()`, and `plot.dm_wavelet()`

    wv <- dm_wavelet(
      x = dm_hourly,
      TreeNum = 1:2,
      source = "raw",
      make_pval = TRUE,
      verbose = FALSE
    )

    wv_summary <- summary(wv)   # summary.dm_wavelet()
    print(wv_summary)           # print.summary.dm_wavelet()

    plot(wv, type = "series")
    plot(wv, type = "average")
    plot(wv, type = "power")
    plot(wv, series = names(wv$results)[1], type = "power")

### `dm_wavelet_reconstruct()` and its methods

    rec_extract <- dm_wavelet_reconstruct(
      wv,
      mode = "extract",
      lower_hours = 20,
      upper_hours = 28,
      only_sig = TRUE
    )

    rec_remove <- dm_wavelet_reconstruct(
      wv,
      mode = "remove",
      lower_hours = 20,
      upper_hours = 28,
      only_sig = TRUE
    )

    rec_extract_summary <- summary(rec_extract)   # summary.dm_wavelet_reconstruct()
    print(rec_extract_summary)                    # print.summary.dm_wavelet_reconstruct()

    plot(rec_extract)
    plot(rec_remove)

### `dm_wavelet_coherence()` and `plot.dm_wavelet_coherence()`

    wc <- dm_wavelet_coherence(
      dm = dm_hourly,
      clim = clim_hourly,
      TreeNum = 1,
      clim_vars = c("temp", "VPD"),
      source = "raw",
      dm_fun = "mean",
      clim_funs = c(temp = "mean", VPD = "mean"),
      loess_span = 0,
      make.pval = TRUE,
      n.sim = 10,
      verbose = FALSE
    )

    plot(wc)
    plot(wc, type = "cross_power")
    plot(wc, type = "average_coherence")
    plot(wc, type = "phase")
    plot(wc, type = "series")

## Closing note

This vignette is intended to be both a workflow guide and a compact reference for the main user-facing functions in the current dendRoAnalyst analysis pipeline. For classroom teaching, the companion demonstration script can be run section by section; for package documentation, this vignette keeps only the lighter sections evaluated and leaves heavy or interactive examples in place for users to run locally.

## For suggestions, comments and questions please contact : sugam.aryal@fau.de