====================================================================================
MerQur - Ziraat_Orman_Su - SCENARIO + RESULT + COMMENTARY (EN, MERGED)
Data: english/datasets/Ziraat_Orman_Su/
Each analysis: SCENARIO + VARIABLE SELECTION, then RESULT (screen) + COMMENTARY.
====================================================================================

#1  Descriptive Statistics
    file: 01_descriptive_forest_inventory.xlsx
  >> SCENARIO (narration):
    In a forest management study we inventoried 300 trees. For each tree we recorded
    species (type), diameter at breast height (dbh_cm), height (height_m), age
    (age_year), biomass (biomass_kg) and a health score (health_score). Before testing
    any hypothesis, we want the overall picture of the stand: what is the mean diameter
    and height, how is biomass distributed, in what range do health scores fall? That is
    why we begin with descriptive statistics; the mean, standard deviation, minimum-maximum
    and distribution summarize the character of the site at a glance.
  >> VARIABLE SELECTION:
    - Variables: dbh_cm
    - Variables: height_m
    - Variables: age_year
    - Variables: biomass_kg
    - Variables: health_score
    - Group/category: type

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    n = 300 trees  |  5 species (beech 67, black pine 60, scots pine 59, calabrian pine 58, oak 56)
    dbh_cm        : mean = 24.59  sd = 14.24  min-max = 5.0 - 80.0
    height_m      : mean = 15.25  sd = 4.84   min-max = 3.0 - 29.5
    age_year      : mean = 66.43  sd = 31.05  min-max = 10 - 119
    biomass_kg    : mean = 164.21 sd = 116.65 min-max = 10 - 600
    health_score  : mean = 3.27   sd = 1.02   min-max = 1 - 5

>> COMMENTARY (narration):
    We start by drawing the overall picture of the stand. Across 300 trees the mean diameter is about
    24.6 cm and the mean height around 15 meters. But what really stands out is the variability: diameter
    ranges from 5 to 80 cm and biomass from 10 to 600 kilograms. Such wide standard deviations tell us this
    is not a uniform plot but a heterogeneous stand holding different age and size classes. The mean age is
    66 years -- a mature forest -- and the health score sits at 3.3 out of 5, moderate-to-good. We have not
    tested any hypothesis yet; but this descriptive table sets the stage for everything that follows -- the
    wide spread in biomass, in particular, already hints at why the normality test will matter.

====================================================================================

#2  Normality Tests
    file: 02_normality_dbh_height_biomass.xlsx
  >> SCENARIO (narration):
    In a sample of 200 trees we examine whether diameter (dbh_cm), height (height_m) and
    biomass (biomass_kg) are normally distributed. The parametric tests we plan to use next
    -- t-tests and ANOVA -- assume that continuous variables follow a normal distribution.
    Shapiro-Wilk and Kolmogorov-Smirnov tests show numerically whether this assumption holds;
    this matters especially for biomass, which in forestry is often right-skewed, and guides us
    to the correct analysis path.
  >> VARIABLE SELECTION:
    - Variables to test: dbh_cm
    - Variables to test: height_m
    - Variables to test: biomass_kg

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    dbh_cm      : Shapiro-Wilk stat = 0.9954  p = 0.804   Normal
    height_m    : Shapiro-Wilk stat = 0.9361  p < 0.001   Not Normal
    biomass_kg  : Shapiro-Wilk stat = 0.8596  p < 0.001   Not Normal

>> COMMENTARY (narration):
    Before moving to parametric tests, we checked whether three variables are normally distributed. The
    result is discriminating: stem diameter (dbh) is almost perfectly normal, with a p-value of 0.80 that
    comfortably supports normality. Height and especially biomass, by contrast, depart significantly from
    normality with p-values near zero. This is no surprise: biomass is typically right-skewed in forests --
    many small trees, a few very large ones. The practical conclusion: we can confidently use parametric
    tests like the t-test and ANOVA on diameter, but for height and biomass it is wiser to consider
    non-parametric alternatives such as Mann-Whitney and Kruskal-Wallis, or transformations / robust methods.

====================================================================================

#3  One-Sample t-Test
    file: 03_one_sample_t_dbh_reference30.xlsx
  >> SCENARIO (narration):
    Here we test whether the mean stem diameter (dbh_cm) of 100 trees in a given stand differs
    significantly from a target reference of 30 cm at the end of the rotation. We have a single
    group and a known comparison value -- no two-group comparison is involved. The one-sample
    t-test is therefore the right choice: we evaluate against a single benchmark whether the
    stand has reached management maturity.
  >> VARIABLE SELECTION:
    - Test variable: dbh_cm
    - Test value (reference): 30

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left): one-sided is more powerful when the direction is known beforehand.
    - Hedges g: small-sample bias-corrected Cohen's d.
    - Effect-size CI: confidence interval around d.
    - Shapiro-Wilk / K-S: normality assumption checks.
    - Descriptives: mean/SD/SE/median/min/max/skewness/kurtosis.
    - Bootstrap CI: distribution-free CI for the mean by resampling.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    t(99) = -3.7606   p < .001 ***   Cohen d = -0.376 (Small)
    Mean dbh = 28.145 cm   (test mu = 30)   95% CI = [27.17, 29.12]
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We compared the mean stem diameter of 100 trees in this stand against a 30 cm reference -- the target for
    rotation maturity. The result is clear: t(99) = -3.76, p below one in a thousand. So the mean diameter --
    28.1 cm -- is significantly lower than the 30 cm target. The confidence interval lies entirely below 30,
    between 27.2 and 29.1. The effect size, Cohen's d = -0.38, is small-to-moderate. In practice this means
    the stand has not yet reached the target diameter maturity; the roughly 2 cm gap may look small but it is
    statistically real, so a bit more time is needed before the rotation is complete.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left): one-sided is more powerful when the direction is known beforehand.
    - Hedges g: small-sample bias-corrected Cohen's d.
    - Effect-size CI: confidence interval around d.
    - Shapiro-Wilk / K-S: normality assumption checks.
    - Descriptives: mean/SD/SE/median/min/max/skewness/kurtosis.
    - Bootstrap CI: distribution-free CI for the mean by resampling.

====================================================================================

#4  Independent Samples t-Test
    file: 04_independent_t_north_south_height.xlsx
  >> SCENARIO (narration):
    We are curious about the effect of slope aspect on tree height. We compare the heights
    (height_m) of trees growing on north-facing versus south-facing aspects. With two
    independent groups and a continuous outcome, the independent-samples t-test is appropriate:
    does the cooler, moister microclimate of the northern aspect significantly change height
    growth compared with the south?
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: height_m
    - Grouping (2 categories): aspect (north, south)

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left).
    - Variance assumption: Student (equal var) / Welch (unequal var — safer) / Auto (Levene decides).
    - Effect sizes: Hedges g, Glass's delta, CLES = P(X>Y).
    - Effect-size CI; per-group Shapiro; Levene & Bartlett homogeneity.
    - Per-group descriptives; Bootstrap CI for the mean difference.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    t(108) = 2.524   p = 0.013 *   Cohen d = 0.483 (Small)
    north: n=50, mean = 17.49 m   |   south: n=60, mean = 16.05 m   |   diff 95% CI = [0.31, 2.57]
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We were curious about the effect of slope aspect on tree height. The results support our hypothesis:
    north-facing trees average 17.5 meters, south-facing ones 16 meters. The roughly one-and-a-half meter
    difference is statistically significant: t(108) = 2.52, p = 0.013. With Cohen's d = 0.48 the effect is
    small-to-moderate. In forestry terms this makes sense: northern aspects offer a cooler, moister
    microclimate with less water stress, which supports height growth. The confidence interval runs from 0.31
    to 2.57 meters -- entirely away from zero. In short, aspect is a genuine driver of height growth on this site.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left).
    - Variance assumption: Student (equal var) / Welch (unequal var — safer) / Auto (Levene decides).
    - Effect sizes: Hedges g, Glass's delta, CLES = P(X>Y).
    - Effect-size CI; per-group Shapiro; Levene & Bartlett homogeneity.
    - Per-group descriptives; Bootstrap CI for the mean difference.

====================================================================================

#5  Paired Samples t-Test
    file: 05_paired_t_dbh_gubreleme.xlsx
  >> SCENARIO (narration):
    In a fertilization trial we track diameter growth on the same trees. For each tree we measured
    diameter before fertilization (dbh_before_cm) and after one growing season (dbh_post_cm). Since
    the measurements come from the same individuals at two time points, the groups are dependent.
    The paired t-test is the correct choice: did fertilization raise per-tree diameter growth by a
    statistically significant amount?
  >> VARIABLE SELECTION:
    - Measurement 1 (before): dbh_before_cm
    - Measurement 2 (after): dbh_post_cm

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left).
    - Effect sizes: Hedges g; d_av (standardized by the average SD).
    - Effect-size CI; pairwise correlation between the two measures.
    - Shapiro / K-S on the differences; descriptives; Bootstrap CI of the mean difference.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Paired t(59) = -19.5399   p < .001 ***   Cohen d_z = -2.523 (Large)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    Here we compared the same trees before and after fertilization -- each tree is its own control. The result
    is striking: t(59) = -19.54, p far below one in a thousand. The negative sign shows the later measurement
    is larger than the earlier one, i.e. diameters increased. What is truly impressive is the effect size:
    Cohen's d_z = -2.52, more than three times the "large" threshold of 0.8 -- the effect of fertilization on
    diameter growth is not just significant, it is enormous. This is the power of a paired design: because we
    compare the same individuals before and after, we eliminate between-tree variability and see the pure
    effect of fertilization very clearly. In practice, the treatment strongly accelerated diameter growth.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction (two-sided / right / left).
    - Effect sizes: Hedges g; d_av (standardized by the average SD).
    - Effect-size CI; pairwise correlation between the two measures.
    - Shapiro / K-S on the differences; descriptives; Bootstrap CI of the mean difference.

====================================================================================

#6  One-Way ANOVA
    file: 06_one_way_anova_type_biomass.xlsx
  >> SCENARIO (narration):
    We investigate whether different tree species (type) differ in biomass (biomass_kg) production.
    Because we need to compare more than two group means at once, one-way ANOVA is appropriate. If
    the F value is significant, at least one species differs from the others; post-hoc (Tukey) tests
    then tell us which species pairs are responsible.
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: biomass_kg
    - Factor (categorical group): type

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - ANOVA variant: Classic (Fisher) or Welch (robust to unequal variances).
    - Effect sizes: omega-squared and epsilon-squared (less biased than eta-squared).
    - Assumptions: Levene, Bartlett, per-group Shapiro.
    - Descriptives per group; post-hoc (Tukey/Duncan/Bonferroni/Scheffe/Games-Howell).

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    F(3, 156) = 33.7295   p < .001 ***   η² = 0.3934
    DECISION: H0 REJECTED (4 species compared)

>> COMMENTARY (narration):
    We compared the biomass production of different tree species. The ANOVA verdict is clear: F(3,156) = 33.73,
    p below one in a thousand. So at least one significant difference exists among the species in biomass. The
    effect size, eta-squared = 0.39, is large; species alone explains 39 percent of the total variability in
    biomass. This is an expected result in forestry: species choice is one of the most fundamental decisions
    determining a stand's production capacity. ANOVA tells us "at least one difference exists" but not between
    which species; the next step would be a post-hoc test such as Tukey to identify which species pairs differ.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - ANOVA variant: Classic (Fisher) or Welch (robust to unequal variances).
    - Effect sizes: omega-squared and epsilon-squared (less biased than eta-squared).
    - Assumptions: Levene, Bartlett, per-group Shapiro.
    - Descriptives per group; post-hoc (Tukey/Duncan/Bonferroni/Scheffe/Games-Howell).

====================================================================================

#7  Two-Way ANOVA
    file: 07_two_way_anova_type_medium.xlsx
  >> SCENARIO (narration):
    In a nursery trial we examine two factors at once: tree species (type) and rearing medium
    (rearing_medium). We want to see how biomass (biomass_kg) is affected by species, by medium, and
    by their interaction. With two categorical factors and one continuous outcome, two-way ANOVA is
    appropriate: the "species x medium" interaction in particular reveals whether a given species
    performs better or worse than expected in a particular medium.
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: biomass_kg
    - Factor 1 (categorical): type
    - Factor 2 (categorical): rearing_medium

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sum-of-squares type (I/II/III): Type III for unbalanced designs with interaction (SPSS default).
    - Post-hoc (Tukey/Bonferroni/Games-Howell) for 3+ level factors.
    - Effect sizes: partial eta-squared, eta-squared, omega-squared.
    - Levene & residual Shapiro; cell and marginal means tables.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    type (main effect)            : F(1,114) = 136.23   p < .001 ***   η²p = 0.54
    rearing_medium (main effect)  : F(2,114) = 155.63   p < .001 ***   η²p = 0.73
    type x rearing_medium         : F(2,114) =   0.49   p = 0.615 ns   η²p = 0.01

>> COMMENTARY (narration):
    In this nursery trial we examined two factors at once: tree species and rearing medium. Both main effects
    came out extremely strong. The rearing medium explains 73 percent of the variability in biomass
    (eta-squared-p = 0.73); species 54 percent. So both which species you choose and the medium you grow it in
    are decisive for your seedling success. But the most instructive part is the interaction: the species x
    medium interaction is non-significant (p = 0.61). This means the best medium is the best medium regardless
    of species; the medium's advantage does not change from species to species -- the effects simply "add up".
    In practice this simplifies nursery management: you pick the most productive medium and that choice holds
    for every species.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sum-of-squares type (I/II/III): Type III for unbalanced designs with interaction (SPSS default).
    - Post-hoc (Tukey/Bonferroni/Games-Howell) for 3+ level factors.
    - Effect sizes: partial eta-squared, eta-squared, omega-squared.
    - Levene & residual Shapiro; cell and marginal means tables.

====================================================================================

#8  Repeated-Measures ANOVA
    file: 08_repeated_anova_yillara_gore_dbh.xlsx
  >> SCENARIO (narration):
    We followed stem diameter on the same trees over four consecutive years: dbh_t0, dbh_t1, dbh_t2
    and dbh_t3. Because the measurements are taken repeatedly from the same individuals, the groups
    are dependent and an independent test would be biased. Repeated-measures ANOVA is therefore
    appropriate: does diameter show a statistically significant growth trend across the years, and
    how does that increase unfold over time?
  >> VARIABLE SELECTION:
    - Repeated measures: dbh_t0
    - Repeated measures: dbh_t1
    - Repeated measures: dbh_t2
    - Repeated measures: dbh_t3

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sphericity correction: Greenhouse-Geisser when Mauchly's test is violated.
    - Mauchly's sphericity test (W, p).
    - Generalized eta-squared (ges) effect size.
    - Post-hoc pairwise (Bonferroni/Holm); descriptives per level.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    F(3, 147) = 265.9963   p < .001 ***   η²p = 0.8444   (n = 50 trees)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We followed stem diameter on the same 50 trees over four consecutive years. The result is the statistical
    fingerprint of growth: F(3,147) = 266, p far below one in a thousand. Partial eta-squared is 0.84 -- an
    extraordinarily large effect; "time", i.e. the years, explains 84 percent of the variability in diameter
    measurements. This is expected but important: because repeated-measures ANOVA tracks the same individuals
    across years, it captures the growth signal very clearly by holding between-tree differences aside. In
    practice this confirms that the stand shows healthy, steady diameter increase -- an active, productive
    growth phase.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sphericity correction: Greenhouse-Geisser when Mauchly's test is violated.
    - Mauchly's sphericity test (W, p).
    - Generalized eta-squared (ges) effect size.
    - Post-hoc pairwise (Bonferroni/Holm); descriptives per level.

====================================================================================

#9  MANOVA
    file: 09_manova_type_height_dbh_biomass.xlsx
  >> SCENARIO (narration):
    We want to examine the effect of tree species (type) not on a single trait but on the whole set
    of height (height_m), diameter (dbh_cm) and biomass (biomass_kg) together. Because these three
    dependent variables are correlated, multivariate analysis of variance (MANOVA) is appropriate
    rather than separate ANOVAs: does species significantly differentiate the combined morphological
    profile? MANOVA also prevents Type-I error inflation.
  >> VARIABLE SELECTION:
    - Dependent variables: height_m
    - Dependent variables: dbh_cm
    - Dependent variables: biomass_kg
    - Factor (categorical group): type

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Reference test for the overall decision: Wilks / Pillai (most robust) / Hotelling-Lawley / Roy.
    - Box's M: equality of covariance matrices across groups.
    - Univariate follow-up ANOVAs (one per dependent variable).
    - Multivariate partial eta-squared; per-group descriptive means.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Wilks' Lambda  ->  F = 131.8669   p < .001 ***
    Dependent variables: height_m, dbh_cm, biomass_kg   |   Factor: type
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    This time we examined the effect of species not on a single trait but on the combined profile of height,
    diameter and biomass. The MANOVA result is very strong: F = 131.87 via Wilks' Lambda, p below one in a
    thousand. So tree species significantly differentiates the joint profile of these three morphological
    traits. The reason we chose MANOVA matters: because these three variables are correlated, running three
    separate ANOVAs would inflate the Type-I error and miss the shared structure among the variables. With a
    single multivariate test we both avoided that risk and reached the conclusion -- "species determines the
    holistic morphology of the tree" -- on solid ground. For detail, one can drill down to each dependent
    variable's individual ANOVA.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Reference test for the overall decision: Wilks / Pillai (most robust) / Hotelling-Lawley / Roy.
    - Box's M: equality of covariance matrices across groups.
    - Univariate follow-up ANOVAs (one per dependent variable).
    - Multivariate partial eta-squared; per-group descriptive means.

====================================================================================

#10  ANCOVA
    file: 10_ancova_treatment_dbh.xlsx
  >> SCENARIO (narration):
    We test the effect of a growth-regulator treatment (treatment) on tree diameter; however, the
    trees already differed in their starting diameter (baseline_dbh_cm). If we compared the post
    diameter (result_dbh_cm) directly, the initial difference would bias the result. We therefore use
    ANCOVA: we enter the baseline diameter as a covariate and hold its effect statistically constant,
    so we see the true, isolated effect of the treatment.
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: result_dbh_cm
    - Covariate (control): baseline_dbh_cm
    - Factor (categorical group): treatment

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sum-of-squares type (I/II/III).
    - Homogeneity-of-regression-slopes test (factor x covariate interaction — the key ANCOVA assumption).
    - Effect sizes: omega-squared, epsilon-squared.
    - Levene & residual Shapiro; Bonferroni post-hoc on adjusted (estimated marginal) means.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    After adjusting for the covariate (baseline diameter), the group difference is SIGNIFICANT.
    Effect size η²p = 0.722 (large)   |   adjusted means reported.
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    Here we tested the effect of a growth-regulator treatment on diameter; but there was a trap: the trees
    already differed in their baseline diameter. If we compared the post diameter directly, that initial
    difference would bias the result. ANCOVA solves exactly this: we entered the baseline diameter as a
    covariate and held its effect statistically constant. The result: even after adjusting for the covariate,
    a significant difference remains between the groups, with a very large effect size of partial eta-squared
    0.72. So the difference we observe is not a by-product of initial differences; it is the true, isolated
    effect of the treatment. ANCOVA is the correct way to fairly compare groups that did not start equal, and
    here it clearly shows the treatment made a strong contribution to diameter growth.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Sum-of-squares type (I/II/III).
    - Homogeneity-of-regression-slopes test (factor x covariate interaction — the key ANCOVA assumption).
    - Effect sizes: omega-squared, epsilon-squared.
    - Levene & residual Shapiro; Bonferroni post-hoc on adjusted (estimated marginal) means.

====================================================================================

#11  Bootstrap Confidence Interval
    file: 11_bootstrap_ci_biomass_median.xlsx
  >> SCENARIO (narration):
    We want to estimate the median of stand biomass (biomass_kg); but because biomass is right-skewed,
    a classical (normality-assuming) confidence interval is not reliable. We therefore use the bootstrap:
    we resample the data thousands of times and build the empirical distribution of the median. This gives
    a robust 95% confidence interval for the biomass median without any distributional assumption.
  >> VARIABLE SELECTION:
    - Variable: biomass_kg
    - Statistic: median

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Column: biomass_kg | Statistic: Median | n = 35 | B = 10,000 replications
    Observed median = 95.0 kg | Bootstrap 95% CI = (75.0, 143.0)

>> COMMENTARY (narration):
    We wanted to estimate the median of stand biomass. Because biomass is right-skewed, a classical
    normality-based confidence interval would not be reliable; so we used the bootstrap. By resampling the data
    ten thousand times we built the empirical distribution of the median. The result: an observed median of 95
    kilograms, and -- most valuable -- a 95% confidence interval saying the true median lies between 75 and 143
    kilograms. The width of this interval is a natural reflection of the high variability in biomass. The beauty
    of the bootstrap is exactly this: with no distributional assumption, we obtained a robust uncertainty interval
    from the data itself.

====================================================================================

#12  Permutation Test
    file: 12_permutation_test_disease_score.xlsx
  >> SCENARIO (narration):
    We test whether two site groups (group) differ in disease score (disease_score) under a small sample
    and uncertain distribution. The permutation test reshuffles the group labels thousands of times and
    computes directly how far the observed difference is from chance. Because it requires no parametric
    assumption, it is the most robust choice in this setting.
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: disease_score
    - Grouping (2 categories): group

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Groups: control (n=28), case (n=25) | disease_score
    DECISION: H0 ACCEPTED (p >= 0.05) -- no significant difference

>> COMMENTARY (narration):
    We tested whether two site groups differ in disease score, under a small sample and uncertain distribution.
    The permutation test is ideal for exactly this situation: it reshuffles the group labels thousands of times
    and computes how ordinary the observed difference is in that random world. The result: the p-value exceeds
    0.05, so the observed difference cannot be distinguished from chance. The statistical decision is clear:
    there is no significant difference in disease score between these two groups. A non-significant result is also
    information -- here, the groups are similar in disease burden.

====================================================================================

#13  Multiple Comparisons (p-value adjustment)
    file: 13_multiple_comparison_conservation_application.xlsx
  >> SCENARIO (narration):
    When we test the effect of several conservation treatments (treatment) on drying percentage (drying_pct)
    separately, performing many comparisons inflates the false-positive risk. Multiple-comparison correction
    (Bonferroni / FDR) adjusts the resulting p-values for the number of comparisons and keeps the Type-I error
    under control, so that "significant" findings are genuinely trustworthy.
  >> VARIABLE SELECTION:
    - Grouping (treatments): treatment
    - Outcome variable: drying_pct
    - Correction: Bonferroni or FDR (Benjamini-Hochberg)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    8 comparisons (each treatment Ui vs control) -> raw p-values before correction
    Bonferroni / Sidak / Holm / BH-FDR / BY-FDR -> TOTAL SIGNIFICANT: 5/8 (all methods)

>> COMMENTARY (narration):
    We compared each of eight conservation treatments against the control separately; but running eight tests
    accumulates false-positive risk. Multiple-comparison correction prevents exactly this. Interestingly, all
    five correction methods -- from the strictest Bonferroni to the modern FDR -- gave the same answer: five of
    the eight treatments still differ significantly from the control even after correction. This shows the result
    is robust; the choice of method did not change it. The practical message is clear: these five treatments
    genuinely work, while three cannot be distinguished from the control. Without correction, we risked reporting
    extra "significant" findings by chance.

====================================================================================

#14  Mann-Whitney U Test
    file: 14_mann_whitney_conservation_area.xlsx
  >> SCENARIO (narration):
    We compare the health scores (health_score) of trees growing in protected versus unprotected areas (area).
    Because the health score is an ordinal/skewed variable and the normality assumption is not met, instead of
    the independent-samples t-test we use its non-parametric counterpart, the Mann-Whitney U test, which
    compares the rank distributions of the two groups.
  >> VARIABLE SELECTION:
    - Dependent (continuous/ordinal) variable: health_score
    - Grouping (2 categories): area

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction; continuity correction.
    - Computation method: auto / exact (precise for small n) / asymptotic.
    - Effect sizes: CLES and Z/sqrt(N) (rank-biserial r already shown).
    - Descriptives (median etc.).

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    U = 1662.50   p = 0.003 **   r = -0.33   (two areas, health_score)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We compared the health scores of trees growing in two different areas. Because the health score is an ordinal,
    skewed variable, we used the non-parametric counterpart of the t-test, the Mann-Whitney U test. The result is
    significant: U = 1662.5, p = 0.003. The effect size r = -0.33 indicates a moderate difference. The two areas
    genuinely differ in tree health. This suggests the conservation status or area conditions have a measurable
    effect on tree health; the rank-based test captured this reliably without needing a normality assumption.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction; continuity correction.
    - Computation method: auto / exact (precise for small n) / asymptotic.
    - Effect sizes: CLES and Z/sqrt(N) (rank-biserial r already shown).
    - Descriptives (median etc.).

====================================================================================

#15  Wilcoxon Signed-Rank Test
    file: 15_wilcoxon_before_post_health.xlsx
  >> SCENARIO (narration):
    We compare the health score of the same trees before (health_before) and after (health_post) a maintenance
    intervention. The measurements are paired, but the health score is not normally distributed, so instead of
    the paired t-test we use its non-parametric counterpart, the Wilcoxon signed-rank test: did the intervention
    produce a statistically significant increase in health?
  >> VARIABLE SELECTION:
    - Measurement 1 (before): health_before
    - Measurement 2 (after): health_post

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction; continuity correction.
    - Zero-difference handling: wilcox (drop) / pratt / zsplit.
    - Descriptives for both measures and their difference.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    W = 0.00   p < .001 ***   r = 0.90   n (non-zero) = 23   (health_before vs health_post)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We compared the health score of the same trees before and after a maintenance intervention. Because health
    score is not normally distributed, we used the non-parametric counterpart of the paired t-test, the Wilcoxon
    signed-rank test. The result is striking: the W statistic is exactly zero, p below one in a thousand. A W of
    zero is a very special case: all trees changed in the same direction -- every single one improved in health.
    The effect size r = 0.90 is enormous. The message is clear: the intervention improved tree health
    consistently and strongly -- not a single tree regressed.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Hypothesis direction; continuity correction.
    - Zero-difference handling: wilcox (drop) / pratt / zsplit.
    - Descriptives for both measures and their difference.

====================================================================================

#16  Kruskal-Wallis Test
    file: 16_kruskal_wallis_soil_development.xlsx
  >> SCENARIO (narration):
    We examine the effect of different soil types (soil_type) on annual seedling development
    (annual_development_cm). There are more than two groups, but the development value is not normally
    distributed; so instead of one-way ANOVA we use its non-parametric counterpart, the Kruskal-Wallis test.
    Is the development ranking of at least one soil type significantly different from the others?
  >> VARIABLE SELECTION:
    - Dependent (continuous/ordinal) variable: annual_development_cm
    - Factor (3+ categories): soil_type

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Epsilon-squared effect size.
    - Dunn post-hoc pairwise comparison (tie-corrected, Bonferroni/Holm).
    - Descriptives per group.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    H(2) = 3.16   p = 0.206 ns   η²_H = 0.03   (3 soil types: sandy/loamy/clayey)
    DECISION: H0 NOT REJECTED

>> COMMENTARY (narration):
    We examined the effect of three soil types on annual seedling development. Because the development value is
    not normally distributed, we used Kruskal-Wallis instead of ANOVA. The result: H(2) = 3.16, p = 0.21 -- not
    significant. The effect size is also very small (0.03). The statistical decision: there is no significant
    difference in seedling development among the three soil types. This too is a valuable finding: at least in
    this trial, soil type did not turn out to be decisive for development; other factors (irrigation, light,
    genetics) may matter more. The non-significance clarifies the "no difference" information.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Epsilon-squared effect size.
    - Dunn post-hoc pairwise comparison (tie-corrected, Bonferroni/Holm).
    - Descriptives per group.

====================================================================================

#17  Friedman Test
    file: 17_friedman_panelist_type_aesthetic.xlsx
  >> SCENARIO (narration):
    In a landscape evaluation, each panelist (panelist_id) rated four tree species -- acacia, poplar, oak and
    pine -- for aesthetics. Because the same evaluators rated all species, the measurements are dependent and
    ordinal. We therefore use the non-parametric counterpart of repeated-measures ANOVA, the Friedman test: is
    the aesthetic-preference ranking significantly different across species?
  >> VARIABLE SELECTION:
    - Repeated measures: acacia
    - Repeated measures: poplar
    - Repeated measures: oak
    - Repeated measures: pine

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Pairwise Wilcoxon signed-rank post-hoc (Bonferroni/Holm).
    - Descriptives per condition.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Friedman χ²(3) = 28.96   p < .001 ***   Kendall W = 0.80   n = 12 panelists
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We examined how twelve panelists rated four tree species -- acacia, poplar, oak, pine -- for aesthetics.
    Because the same evaluators rated all species, we used the non-parametric counterpart of repeated-measures
    ANOVA, the Friedman test. The result is highly significant: chi-square 28.96, p below one in a thousand.
    What is most striking is Kendall's W = 0.80: this indicates very high agreement among the panelists. So the
    evaluators largely agree on aesthetic preference -- certain species are consistently preferred over others.
    In landscape and urban forestry this matters: species choice, though it seems subjective, actually rests on
    a strong shared aesthetic judgment.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Pairwise Wilcoxon signed-rank post-hoc (Bonferroni/Holm).
    - Descriptives per condition.

====================================================================================

#18  Binomial Test
    file: 18_binomial_seedling_survival_ratio075.xlsx
  >> SCENARIO (narration):
    On an afforestation site we test the survival rate of the seedlings we planted. Whether each seedling took
    (experienced: survived/died) was recorded, and the target survival rate is set at 0.75. With a single binary
    outcome and a known reference proportion, the binomial test is appropriate: is the observed survival rate
    significantly different from the targeted 0.75?
  >> VARIABLE SELECTION:
    - Binary outcome variable: experienced (survived/died)
    - Test proportion: 0.75

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Observed proportion = 0.725 | Expected (target) = 0.750 | p = 0.073 ns
    DECISION: H0 ACCEPTED -- proportion not significantly different from 0.75

>> COMMENTARY (narration):
    We compared the survival rate of seedlings on an afforestation site against the target of 0.75. The observed
    rate is 0.725, i.e. 72.5 percent. The binomial test asks: is this below target, or just chance fluctuation?
    The answer, with p = 0.073, is clear: the difference is not statistically significant. So the observed
    survival rate is consistent with the targeted 75 percent; the small gap can be explained by chance. In
    practice this is good news: the afforestation can be considered to have met its survival target. Note the
    p-value of 0.073 is close to the threshold -- with a larger sample the difference could become significant,
    so it is sensible to keep monitoring.

====================================================================================

#19  Sign Test
    file: 19_sign_test_pruning_crown_diameter.xlsx
  >> SCENARIO (narration):
    We compare the crown diameter of the same trees before (crown_diameter_before_m) and after
    (crown_diameter_post_m) pruning. The data are paired, but the distribution of differences is not symmetric;
    in this case, which even strains the symmetry assumption of Wilcoxon, the sign test -- which considers only
    the direction of increase/decrease -- is the most robust choice: did pruning systematically change crown
    diameter?
  >> VARIABLE SELECTION:
    - Measurement 1 (before): crown_diameter_before_m
    - Measurement 2 (after): crown_diameter_post_m

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Positive (before > after): 28 | Negative (before < after): 2 | p < .001 ***
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We compared the crown diameter of the same trees before and after pruning. Because the distribution of
    differences is not symmetric, we used the sign test -- which looks only at the direction of change -- instead
    of Wilcoxon. The table is very clear: in 28 of 30 trees the crown diameter shrank after pruning, and increased
    in only 2. The p-value is below one in a thousand. So pruning systematically and significantly reduced crown
    diameter -- which is, of course, the purpose of pruning. The power of the sign test is here: with no
    distributional assumption, relying only on counting "decreased/increased", we reached a very safe conclusion.

====================================================================================

#20  Runs Test
    file: 20_runs_test_fire_direction.xlsx
  >> SCENARIO (narration):
    In a forest-fire monitoring record we wonder whether the daily spread direction (direction_binary: two
    directions such as north/south) follows a random pattern or a clustered one. To test the randomness of the
    sequence, the runs test is appropriate: are consecutive same-direction runs fewer or more than expected --
    that is, is there a directional tendency or autocorrelation in the spread?
  >> VARIABLE SELECTION:
    - Binary sequence variable: direction_binary

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Number of runs BELOW expected | p < .001 *** -> CLUSTERING detected
    DECISION: H0 REJECTED (sequence is not random)

>> COMMENTARY (narration):
    We tested whether the daily spread direction of a forest fire follows a random pattern or a structured one.
    The runs test looks at the number of consecutive same-direction runs. The result: the number of runs is
    clearly below expected, p below one in a thousand. This means clustering: the same direction repeats in a
    row, so the fire direction is not random from day to day but shows a persistent tendency. This is a very
    important practical finding: if spread carries directional autocorrelation, fire behavior becomes
    predictable -- directly usable information for intervention and resource-allocation planning.

====================================================================================

#21  Chi-Square Independence Test
    file: 21_chisquare_independence_type_disease.xlsx
  >> SCENARIO (narration):
    We examine whether there is a relationship between tree species (type) and disease occurrence
    (disease_present: yes/no). To test the independence of two categorical variables, the chi-square test of
    independence is appropriate: are certain species more prone to disease than others, or are species and
    disease independent of each other?
  >> VARIABLE SELECTION:
    - Row variable (categorical): type
    - Column variable (categorical): disease_present

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Yates continuity correction (2x2 tables).
    - G-test (likelihood-ratio chi-square) alternative.
    - Effect sizes: phi (2x2) and contingency coefficient C (besides Cramer's V).
    - Expected-counts table; standardized residuals (|>2| flags the deviating cell).

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    χ²(3) = 17.42   p < .001 ***   Cramer's V = 0.21 (Medium)   (type x disease_present)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We examined whether there is a relationship between tree species and disease occurrence with the chi-square
    test of independence. The result is significant: chi-square 17.42, p below one in a thousand. So species and
    disease are not independent; some species are more prone to disease than others. Cramer's V = 0.21 means the
    relationship is moderate -- not weak, but not all-determining either. In forestry this matters: if disease
    risk depends on species choice, management decisions that favor resistant species can reduce the disease burden.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Yates continuity correction (2x2 tables).
    - G-test (likelihood-ratio chi-square) alternative.
    - Effect sizes: phi (2x2) and contingency coefficient C (besides Cramer's V).
    - Expected-counts table; standardized residuals (|>2| flags the deviating cell).

====================================================================================

#22  Chi-Square Goodness-of-Fit
    file: 22_chisquare_goodnessfit_type_equal_distribution.xlsx
  >> SCENARIO (narration):
    We test whether tree species (type) are distributed equally in a natural stand, or whether certain species
    dominate. To test the observed distribution of a single categorical variable against an expected (equal)
    distribution, the chi-square goodness-of-fit test is appropriate: is the species distribution homogeneous,
    or is there a significant departure?
  >> VARIABLE SELECTION:
    - Categorical variable: type
    - Expected distribution: equal (uniform)

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Effect sizes: Cohen's w and Cramer's V.
    - G-test (likelihood ratio) alternative.
    - Standardized residuals per category (|>2| = notable deviation).

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    χ²(3) = 5.00   p = 0.172 ns   (n = 200, k = 4 species, expected: equal distribution)
    DECISION: H0 NOT REJECTED

>> COMMENTARY (narration):
    We tested whether four tree species are distributed equally in this natural stand. The chi-square
    goodness-of-fit test compared the observed species distribution against an "all equal" expectation. The
    result: chi-square 5.0, p = 0.17 -- not significant. So there is no statistically significant departure from
    an equal distribution; the stand is quite balanced in species composition. This is a positive picture for
    biodiversity: no single species dominates, and species are represented in similar proportions.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Effect sizes: Cohen's w and Cramer's V.
    - G-test (likelihood ratio) alternative.
    - Standardized residuals per category (|>2| = notable deviation).

====================================================================================

#23  Fisher's Exact Test
    file: 23_fisher_exact_small_sample_seedling.xlsx
  >> SCENARIO (narration):
    In a small seedling trial we examine the relationship between two species (type) and the success outcome
    (result: took/failed). Because the sample is small and some cell frequencies are below 5, the chi-square is
    not reliable; in this case Fisher's exact test, which computes the exact probability directly, is appropriate:
    is there a significant dependence between species and seedling establishment success?
  >> VARIABLE SELECTION:
    - Row variable (categorical): type
    - Column variable (categorical): result

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Fisher exact p = 0.025 *   Odds Ratio = 0.125   (95% CI: [0.024, 0.657])   (type x result)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    In a small seedling trial we examined the relationship between two species and seedling establishment
    success. Because the sample is small and some cells fall below 5, we used Fisher's exact test, which computes
    the exact probability, instead of chi-square. The result is significant: p = 0.025. The odds ratio is 0.125,
    so the chance of seedling establishment differs markedly between species; the confidence interval lies
    entirely below 1 (0.024-0.657), which clarifies the direction of the relationship. The value of Fisher's test
    in small samples is exactly this: instead of an approximate test, we made a safe decision with exact probability.

====================================================================================

#24  McNemar Test
    file: 24_mcnemar_old_new_diagnosis.xlsx
  >> SCENARIO (narration):
    We compare the diagnostic outcomes of an old method (old_method) and a new method (new_method) for detecting
    tree disease on the same trees. Because the measurements are paired and binary (diseased/healthy), the McNemar
    test is appropriate: are the discordant cases (one positive, the other negative) balanced, or does the new
    method systematically diagnose more/less?
  >> VARIABLE SELECTION:
    - Measurement 1 (binary): old_method
    - Measurement 2 (binary): new_method

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Discordant: b = 8, c = 12 | min(b,c) = 8 | p = 0.503 ns   (old_method vs new_method)
    DECISION: H0 NOT REJECTED

>> COMMENTARY (narration):
    We compared the diagnostic outcomes of an old and a new method for detecting tree disease on the same trees.
    The McNemar test looks only at the DISCORDANT cases: there are 20 trees where one method said positive and the
    other negative (8 + 12). If the new method systematically diagnosed more or fewer, these discordances would be
    unbalanced. But with p = 0.50 the difference is not significant: an 8-versus-12 split is indistinguishable from
    chance. Conclusion: the two methods are statistically equivalent in diagnosis; there is no systematic shift in
    the new method. This suggests the new method can safely replace the old one.

====================================================================================

#25  Cohen's Kappa
    file: 25_kappa_two_expert_health.xlsx
  >> SCENARIO (narration):
    We measure the agreement between two experts (expert_A, expert_B) who independently classify the health status
    of the same trees. Simple percent agreement also includes agreement by chance; Cohen's kappa corrects for chance
    and gives the true agreement: how consistent are the experts in their health assessment -- is kappa at a "good"
    level (0.6+)?
  >> VARIABLE SELECTION:
    - Rater 1: expert_A
    - Rater 2: expert_B

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Cohen's κ = 0.777 (substantial)   p_o = 0.85   p_e = 0.33   (95% CI: 0.66+)
    DECISION: H0 REJECTED (agreement beyond chance)

>> COMMENTARY (narration):
    We measured the agreement between two experts independently classifying the health status of the same trees.
    The raw agreement appears to be 85 percent; but part of that could be due to chance. Cohen's kappa corrects
    for chance and gives the true agreement: kappa = 0.78. On the kappa scale this is "substantial" -- the experts
    are largely consistent in their health assessment. In practice this is reassuring: the assessment protocol is
    reliable, and whoever measures gives a similar result. In health-monitoring programs, such high inter-rater
    agreement is a basic requirement for data reliability.

====================================================================================

#26  Cochran-Mantel-Haenszel (CMH)
    file: 26_cmh_conservation_disease_region.xlsx
  >> SCENARIO (narration):
    We examine the relationship between being in a conservation area (conservation_area) and disease (disease),
    controlling for the effect of region (region). If region is a confounder that affects both protection status
    and disease risk, the relationship must be handled separately within region strata. The CMH test combines these
    stratified 2x2 tables and shows whether the conservation-disease association is a common effect independent of
    region.
  >> VARIABLE SELECTION:
    - Row (categorical): conservation_area
    - Column (categorical): disease
    - Stratum/layer (control): region

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    CMH χ²(1) = 20.95   p < .001 ***   Common Odds Ratio (MH) = 0.362   (stratified by region)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We examined the relationship between being in a conservation area and disease, controlling for region. Region
    could be a confounder affecting both protection status and disease risk, so we handled the relationship within
    region strata. The CMH result is highly significant: chi-square 20.95, p below one in a thousand. The common
    odds ratio is 0.36, meaning that EVEN AFTER accounting for regional differences, being in a conservation area
    markedly reduces disease risk (OR below 1 = protective effect). This shows the conservation status has a real,
    region-independent effect -- one of the most valuable findings.

====================================================================================

#27  Log-Linear Models
    file: 27_log_linear_type_region_disease.xlsx
  >> SCENARIO (narration):
    We want to examine the main-effect and interaction structure in the contingency table formed jointly by three
    categorical variables -- tree species (type), region (region) and disease (disease). Chi-square compares only two
    variables; log-linear models, by modeling the cell frequencies of three or more variables via exp(B'X), reveal
    which interactions (e.g. species x disease) are significant.
  >> VARIABLE SELECTION:
    - Categorical variable: type
    - Categorical variable: region
    - Categorical variable: disease

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Log-Linear Model (order = 2)   AIC = 99.66   Pearson χ² = 10.17
    log(cell_count) = β'X   |   variables: type, region, disease

>> COMMENTARY (narration):
    We examined the structure of the contingency table formed jointly by three categorical variables -- species,
    region and disease -- with a log-linear model. Chi-square can only test pairwise relationships; the log-linear
    model, by modeling the cell frequencies of three variables at once, shows which interactions matter. The
    two-way (order-2) model fits the data reasonably (Pearson chi-square 10.17, AIC 99.66); this suggests that
    pairwise dependencies (e.g. species-disease) explain the table and no complex three-way interaction is needed.
    In practice this means the disease pattern can be explained by the simple pairwise effects of species and region.

====================================================================================

#28  Cross-Tabulation Analysis
    file: 28_cross_tablo_age_type_preference.xlsx
  >> SCENARIO (narration):
    We examine the relationship between forest visitors' age group (age_group) and their preferred forest type
    (preference_type) with a cross-tabulation. The cross-tab shows the joint distribution of two categorical
    variables with row/column percentages and tests the significance of the association with chi-square: do younger
    and older visitors prefer different forest types?
  >> VARIABLE SELECTION:
    - Row variable (categorical): age_group
    - Column variable (categorical): preference_type

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    χ²(6) = 22.30   p = 0.001 ***   Cramer's V = 0.15 (Medium)   (age_group x preference_type)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We examined the relationship between forest visitors' age group and their preferred forest type with a
    cross-tabulation. The result is significant: chi-square 22.30, p = 0.001. So age and forest-type preference
    are not independent; different age groups lean toward different forest types. Cramer's V = 0.15 means the
    relationship is moderate-to-weak -- real, but not the sole determinant. This is useful for recreation
    planning: when designing areas that target different age groups, knowing which age group prefers which forest
    type translates directly into application.

====================================================================================

#29  Multiple-Response Frequency
    file: 29_mr_frequency_forest_activity.xlsx
  >> SCENARIO (narration):
    In a recreation survey, participants were asked which activities they do in the forest; each person could check
    multiple options (walking, camping, bird observation, photography, picnic, fungus, cycling). Each activity is a
    separate binary (0/1) column. Multiple-response frequency analysis summarizes these options together: which
    activity is most preferred, and what are the response and case percentages?
  >> VARIABLE SELECTION:
    - Multiple-response (binary) columns: activity_walking, activity_camping, activity_bird_observation,
      activity_photography, activity_picnic, activity_fungus, activity_cycling

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Total cases = 350 | Respondents = 339 (96.9%) | Total responses = 935
    Walking 70.2% (238) > Picnic 53.4% > Photography 47.8% > Cycling 32.5% > Camping 28.9% ...

>> COMMENTARY (narration):
    In a recreation survey we asked participants which activities they do in the forest; everyone could check
    multiple options. Multiple-response frequency analysis summarizes these overlapping options together. The 339
    respondents marked a total of 935 activities -- about three activities per person. Walking is the clear leader:
    70 percent of participants do it. It is followed by picnicking (53%) and photography (48%). Foraging and bird
    watching are more niche activities. This table directly tells us where to prioritize when designing forest
    recreation areas: walking trails and picnic areas yield the highest return.

====================================================================================

#30  Multiple-Response Crosstab
    file: 30_mr_categorical_activity_region.xlsx
  >> SCENARIO (narration):
    This time we cross the multiple-response forest activities (walking, camping, picnic, photography) by region
    (region). A multiple-response x categorical cross-tab shows which activities stand out in each region: e.g. is
    camping dominant in the mountains and walking on the coast? Is there a relationship between region and activity
    preference?
  >> VARIABLE SELECTION:
    - Multiple-response (binary) columns: activity_walking, activity_camping, activity_picnic, activity_photography
    - Categorical breakdown: region

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Total cases = 300 | Group: region (Aladag, Belgrad, Kackar)
    case-percentage table of each region x activity (walking/camping/picnic/photography)

>> COMMENTARY (narration):
    This time we crossed the multiple-response activities by region. The multiple-response x categorical table
    shows which activities stand out in each region. This way we capture regional patterns that a single overall
    average would hide: e.g. camping may dominate in a mountainous region while walking dominates in an urban
    forest. The fact that activity preference varies by region shows that a "one-size-fits-all" recreation plan is
    not enough; each area should be designed for its own visitor profile. This is a practical, actionable finding
    for directing resources to the right place.

====================================================================================

#31  Multiple-Response x Multiple-Response Crosstab
    file: 31_mr_mr_activity_equipment.xlsx
  >> SCENARIO (narration):
    Now we cross two multiple-response sets: the activities participants do (walking, camping, picnic) and the
    equipment they bring (tent, back bag, photography machine, camera, GPS). A multiple-response x multiple-response
    table shows which activity goes with which equipment: e.g. do campers carry tent + GPS together? What is the
    activity-equipment co-occurrence pattern?
  >> VARIABLE SELECTION:
    - Set 1 (multiple-response): activity_walking, activity_camping, activity_picnic
    - Set 2 (multiple-response): crew_tent, crew_back_bag, crew_photography_machine, crew_camera, crew_gps

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Set 1 (rows): activity (walking/camping/picnic) | Set 2 (cols): crew (tent/bag/machine/camera/GPS)
    Total cases = 250 | activity x equipment co-occurrence table

>> COMMENTARY (narration):
    Here we crossed two multiple-response sets: the activities participants do and the equipment they bring. The
    multiple-response x multiple-response table shows which activity goes with which equipment. This is richer than
    a simple frequency: instead of just "camping is popular" and "tents are popular", we can see that campers carry
    tent and GPS together, and photographers carry camera and machine. The practical value is large: if we know the
    activity profile of an area, we can anticipate which equipment visitors will need and plan infrastructure
    (campsites, charging points, rentals) accordingly.

====================================================================================

#32  Cochran's Q Test
    file: 32_cochrans_q_3saha_safety.xlsx
  >> SCENARIO (narration):
    We had the same visitors (visitor_id) evaluate whether three different recreation sites felt safe
    (site_A_safe, site_B_safe, site_C_safe -- each 0/1). Because the same subjects rated all three sites and the
    response is binary, Cochran's Q test is appropriate (the binary counterpart of Friedman): are the "felt-safe"
    rates of the three sites significantly different from one another?
  >> VARIABLE SELECTION:
    - Binary (0/1) columns: site_A_safe, site_B_safe, site_C_safe

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Cochran's Q = 35.51   p < .001 ***   (3 sites: site_A/B/C safe, 0/1)
    DECISION: H0 REJECTED

>> COMMENTARY (narration):
    We had the same visitors evaluate whether three recreation sites felt safe. Because the same subjects rated all
    three sites and the response is binary (yes/no), we used Cochran's Q -- the binary counterpart of Friedman. The
    result is highly significant: Q = 35.51, p below one in a thousand. So the felt-safe rates of the three sites
    genuinely differ; at least one site is perceived as markedly safer (or less safe) than the others. This is a
    management-oriented finding that directly indicates which site needs priority intervention to improve its safety
    perception.

====================================================================================

#33  Correlation Analysis
    file: 33_correlation_forest_measurements.xlsx
  >> SCENARIO (narration):
    We examine the relationships among forest measurement variables: age (age_year), diameter (dbh_cm), height
    (height_m), biomass (biomass_kg) and crown diameter (crown_m). Correlation analysis summarizes the pairwise
    relationships (Pearson r) of these continuous variables and their direction/strength: is there the expected
    strong positive relationship between diameter and biomass, and which measurements best predict one another?
  >> VARIABLE SELECTION:
    - Variables: age_year, dbh_cm, height_m, biomass_kg, crown_m

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - p-values are now produced for ALL methods (Pearson/Spearman/Kendall), not only Pearson.
    - Hypothesis direction (two-sided / right / left).
    - Multiple-comparison p-adjustment across pairs: Bonferroni / Holm / FDR (Benjamini-Hochberg).
    - Confidence interval for r via Fisher z (Pearson/Spearman) or Kendall SE.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    STRONG relationships (|r| >= 0.75) -- 10 pairs:
    age-height r=0.96 | dbh-biomass r=0.956 | height-biomass r=0.954 | dbh-crown r=0.948
    dbh-height r=0.946 | age-biomass r=0.928   (all p < .0001)

>> COMMENTARY (narration):
    We examined the relationships among forest measurement variables and the table reads like a textbook example.
    All pairwise relationships are very strong and positive: age with height 0.96, diameter with biomass 0.956,
    height with biomass 0.954. All are p far below one in a thousand. This is expected but important: in a tree,
    diameter, height, age and biomass grow together -- as the tree grows, all increase. The especially strong
    diameter-biomass relationship is a practical opportunity: biomass is hard to measure (you must fell the tree),
    but diameter is easy. Thanks to this high correlation, we can estimate biomass non-destructively with
    diameter-based allometric equations.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - p-values are now produced for ALL methods (Pearson/Spearman/Kendall), not only Pearson.
    - Hypothesis direction (two-sided / right / left).
    - Multiple-comparison p-adjustment across pairs: Bonferroni / Holm / FDR (Benjamini-Hochberg).
    - Confidence interval for r via Fisher z (Pearson/Spearman) or Kendall SE.

====================================================================================

#34  Bland-Altman Analysis
    file: 34_bland_altman_clinometer_laser.xlsx
  >> SCENARIO (narration):
    We compare two methods of measuring tree height: the classical clinometer (clinometer_m) and a laser hypsometer
    (laser_m). Correlation tells us the two methods are "related" but not whether they "agree". Bland-Altman analysis
    plots the measurement differences against the mean and shows the systematic bias and limits of agreement: can the
    laser method safely replace the clinometer?
  >> VARIABLE SELECTION:
    - Method 1: clinometer_m
    - Method 2: laser_m

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Bias (mean difference) = 0.39 m   95% CI = (0.24, 0.53)   Bias=0 test: t = 5.27, p < .001
    -> systematic difference EXISTS

>> COMMENTARY (narration):
    We compared two methods of measuring tree height -- the clinometer and a laser. Correlation tells us the two
    methods are "related" but not whether they "agree"; that is what Bland-Altman tests. The result is important:
    the mean difference (bias) is 0.39 meters, and this difference is statistically different from zero (t = 5.27,
    p < .001). So there is a systematic offset between the two methods -- one consistently measures about 40 cm
    more/less than the other. Practical conclusion: these two methods cannot be used interchangeably as-is; if
    switching to the laser, this fixed 0.39 m difference must be accounted for as a correction factor. Bland-Altman
    revealed this critical detail that correlation would have hidden.

====================================================================================

#35  Effect Size (Cohen's d)
    file: 35_effect_size_cultivation_method.xlsx
  >> SCENARIO (narration):
    We assess whether the effect of two seedling cultivation methods (method) on seedling height (height_cm) is not
    only "significant" but also "large". The p-value tells us whether a difference exists; the effect size (Cohen's d)
    gives the practical magnitude of the difference in standard-deviation units: d = 0.2 small, 0.5 medium, 0.8 large.
    Is the method difference meaningful in forestry practice?
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: height_cm
    - Grouping (2 categories): method

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Cohen's d = 1.66 (Large)   95% CI [1.08, 2.25]   (height_cm, method: 2 groups)
    DECISION: very large effect

>> COMMENTARY (narration):
    We assessed whether the effect of two seedling cultivation methods on seedling height is not only significant
    but large in practice. The p-value tells whether a difference exists; the effect size gives the true magnitude.
    Cohen's d = 1.66 -- more than twice the "large" threshold of 0.8, with a confidence interval entirely above 1
    (1.08-2.25). So the difference between the two methods is both significant and enormous: one method raises
    seedling height by about 1.7 standard deviations over the other. In forestry practice this shows that method
    choice genuinely matters -- we did not fall into the "statistically significant but practically trivial" trap.

====================================================================================

#36  Canonical Correlation (CCA)
    file: 36_cca_fizyoloji_growth.xlsx
  >> SCENARIO (narration):
    We examine the joint variation of two multivariate measurement sets: a physiological set (photosynthesis, stomata
    count, leaf area, chloroplast density, water potential) and a growth set (annual diameter increase, annual height
    increase, side-branch count, inner-growth index). Canonical correlation (CCA) finds the linear combinations that
    most relate the two sets: how much, and along which axis, does the tree's physiological capacity predict its growth
    performance?
  >> VARIABLE SELECTION:
    - X set (physiology): photosynthesis, stomata_count, leaf_area, chloroplast_density, water_potential
    - Y set (growth): annual_dbh_increase, annual_height_increase, side_branch_count, inner_growth_index

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    CC1: r = 0.90 (r²=0.81)  χ²(20)=463.6  p < .001 ***
    CC2: r = 0.72 (r²=0.52)  χ²(12)=142.1  p < .001 ***   |   CC3, CC4: non-significant
    (X = physiology 5 vars, Y = growth 4 vars)

>> COMMENTARY (narration):
    We examined two multivariate sets together -- the tree's physiological capacity and its growth performance.
    Canonical correlation finds the linear combinations that most relate the two sets. The result is very strong:
    the first canonical dimension, r = 0.90, links physiology and growth almost perfectly; the second dimension is
    also significant (r = 0.72). The third and fourth are noise. So the tree's physiology (photosynthesis, stomata,
    leaf area, water potential) is strongly related to its growth (diameter/height increase, branch count) along
    two independent axes. This shows physiological measurements are valuable for predicting growth potential -- that
    early physiological screening could forecast future growth performance.

====================================================================================

#37  Correspondence Analysis
    file: 37_correspondence_type_region.xlsx
  >> SCENARIO (narration):
    We want to examine the relationship between tree species (type) and region (region) on a visual map. Correspondence
    analysis (CA) positions the contingency table of categorical variables in a two-dimensional space; species-region
    pairs that appear close together tend to co-occur. Which species are associated with which regions -- what is the
    ecological distribution pattern?
  >> VARIABLE SELECTION:
    - Row variable (categorical): type
    - Column variable (categorical): region

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    N = 600 | Total inertia = 0.683 | Dimensions = 2   (type x region)

>> COMMENTARY (narration):
    We examined the relationship between tree species and region on a visual map. Correspondence analysis places
    the contingency table of categorical variables into a two-dimensional space; species-region pairs that appear
    close tend to co-occur. The total inertia 0.683 means there is a notable association between species and region,
    summarized in two dimensions. In practice this map is instructive: it shows at a glance which species is
    associated with which region -- the ecological distribution pattern. These species-region proximities can be
    used directly to plan natural distribution ranges and region-appropriate species selection in afforestation.

====================================================================================

#38  Variable Clustering (VARCLUS)
    file: 38_varclus_forest_measurements_18_item.xlsx
  >> SCENARIO (narration):
    We have 18 trunk/leaf/root measurements (trunk_M1..M6, leaf_M1..M6, root_M1..M6). Most of these variables are
    highly correlated with one another; putting them all in one model causes multicollinearity. Variable clustering
    (VARCLUS) groups the measurements by their correlation structure and selects a representative from each cluster;
    this reduces 18 variables to a few independent dimensions for data reduction.
  >> VARIABLE SELECTION:
    - Numeric variables: trunk_M1..trunk_M6, leaf_M1..leaf_M6, root_M1..root_M6 (18 measurements)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    18 measurements (trunk/leaf/root M1-M6) -> 3 clusters | n = 250
    one representative selected per cluster (e.g. trunk_M6, ...)

>> COMMENTARY (narration):
    We had 18 trunk, leaf and root measurements; most are highly correlated, and putting them all in one model
    would cause multicollinearity. Variable clustering split these 18 measurements into three natural clusters by
    their correlation structure and selected the variable that best represents each cluster. This reduced 18
    variables to three independent dimensions. The practical benefit is clear: in later analyses, instead of using
    all 18 measurements, we can use the three selected representatives -- avoiding multicollinearity and simplifying
    the model, with minimal information loss. This is a classic, powerful way to do data reduction.

====================================================================================

#39  Linear Regression
    file: 39_lineer_regression_biomass.xlsx
  >> SCENARIO (narration):
    We want to build a model that predicts tree biomass (biomass_kg). Candidate predictors: diameter (dbh_cm), height
    (height_m), age (age_year) and soil pH (soil_pH). Multiple linear regression gives the contribution of each
    predictor to biomass (the B coefficient) while holding the others constant, and shows the variance the model
    explains (R2): is biomass driven mostly by diameter, and how strong is the model?
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: biomass_kg
    - Predictors: dbh_cm, height_m, age_year, soil_pH

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Robust standard errors (HC0-HC3): heteroscedasticity-robust SE; HC3 recommended for small n.
    - Standardized (beta) coefficients to compare relative effect.
    - (Diagnostics VIF, Durbin-Watson, Breusch-Pagan, residual Shapiro are already reported.)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    R² = 0.954   Adj. R² = 0.953   F = 1523.4   p < .001 ***   (n = 300)
    biomass_kg ~ dbh_cm + height_m + age_year + soil_pH

>> COMMENTARY (narration):
    We built a model to predict tree biomass. The result is extraordinarily strong: the model explains 95 percent
    of the variability in biomass (R² = 0.954). F = 1523, p far below one in a thousand -- the model as a whole is
    extremely significant. This high R² is expected in forestry: biomass is almost a deterministic function of the
    tree's diameter, height and age. The practical value is very high: with this equation we can estimate a tree's
    biomass to 95 percent accuracy without felling it, just by measuring its diameter, height and age. From carbon
    stock calculation to harvest planning, such allometric models are fundamental tools of forestry.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Robust standard errors (HC0-HC3): heteroscedasticity-robust SE; HC3 recommended for small n.
    - Standardized (beta) coefficients to compare relative effect.
    - (Diagnostics VIF, Durbin-Watson, Breusch-Pagan, residual Shapiro are already reported.)

====================================================================================

#40  Logistic Regression
    file: 40_logistic_regression_disease.xlsx
  >> SCENARIO (narration):
    We look for the factors that predict whether a tree will become diseased (disease: 0/1): age (age), stress score
    (stress_score) and soil moisture (moisture_pct). Because the outcome is binary, logistic regression is appropriate:
    it gives each predictor's effect on the probability as an odds ratio (OR). Does disease probability increase
    significantly as the stress score rises, and is moisture protective?
  >> VARIABLE SELECTION:
    - Dependent (binary) variable: disease
    - Predictors: age, stress_score, moisture_pct

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Pseudo-R squared: Cox-Snell and Nagelkerke (besides McFadden).
    - Classification metrics: accuracy / sensitivity / specificity / AUC (cutoff 0.5).
    - Hosmer-Lemeshow goodness-of-fit test; VIF for predictors.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Pseudo R² = 0.21   model significant (p < .001)   n = 400
    disease ~ age + stress_score + moisture_pct

>> COMMENTARY (narration):
    We used logistic regression to find the factors predicting whether a tree becomes diseased: age, stress score
    and soil moisture. The model is significant overall (p < .001) and the McFadden pseudo-R² is 0.21 -- a good fit
    for logistic models (not to be confused with classical R²; in logistic, 0.2-0.4 is already strong). Among the
    predictors, the stress score in particular stands out as the most prominent factor raising disease probability,
    while moisture is protective. In practice this model can act as an early-warning tool: trees under high stress
    and low moisture can be prioritized for monitoring against disease.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Pseudo-R squared: Cox-Snell and Nagelkerke (besides McFadden).
    - Classification metrics: accuracy / sensitivity / specificity / AUC (cutoff 0.5).
    - Hosmer-Lemeshow goodness-of-fit test; VIF for predictors.

====================================================================================

#41  Count Regression (Poisson / Negative Binomial)
    file: 41_poisson_negbin_insect_sayim.xlsx
  >> SCENARIO (narration):
    We model the insect count (insect_count) in each experimental plot. Because count data consist of non-negative
    integers, linear regression is not appropriate; Poisson / Negative Binomial regression is. The predictors are plant
    richness (plant_richness), temperature (temperature_C) and moisture (moisture_pct). exp(B) = the incidence rate ratio
    (IRR): by how many times does the insect count increase as temperature rises? (If overdispersion is present, use
    Negative Binomial.)
  >> VARIABLE SELECTION:
    - Dependent (count) variable: insect_count
    - Predictors: plant_richness, temperature_C, moisture_pct

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Poisson GLM   AIC = 1392.4   Deviance = 246.9   |   IRR table (exp(coef))
    insect_count ~ plant_richness + temperature_C + moisture_pct

>> COMMENTARY (narration):
    We modeled the insect count in each plot. Because count data are non-negative integers, we used Poisson
    regression rather than linear regression. The model's most intuitive output is each predictor's incidence rate
    ratio (IRR = exp(coef)): by how many times the insect count changes when a variable rises by one unit.
    Temperature stands out as the most prominent factor raising insect density; plant richness and moisture also
    contribute. In practice this helps predict which plots are at risk of insect infestation: warm, suitably moist
    areas gain priority in monitoring. (If overdispersion is present, switching to Negative Binomial is more correct.)

====================================================================================

#42  Multinomial Logistic Regression
    file: 42_multinomial_forest_type_preference.xlsx
  >> SCENARIO (narration):
    We try to predict visitors' preferred forest type (preference_forest -- an unordered 3+ category) from age (age)
    and years of education (education_year). Because the outcome has more than two unordered categories, multinomial
    logistic regression is appropriate: it builds a separate equation for each non-reference category and shows how
    age/education change the preference probability.
  >> VARIABLE SELECTION:
    - Dependent (3+ unordered categories): preference_forest
    - Predictors: age, education_year

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Multinomial Logit   AIC = 600.3   |   reference class: 'leaved'
    preference_forest (3+ unordered) ~ age + education_year

>> COMMENTARY (narration):
    We tried to predict visitors' preferred forest type -- an unordered, 3+ category outcome -- from age and
    education. Because the outcome has more than two unordered categories, we used multinomial logistic regression,
    which builds a separate equation for each non-reference category. With "leaved" as reference, we read how age
    and education change the odds of preferring each other forest type. In practice this enables demographic
    targeting in recreation planning: e.g. if higher education tends toward a particular forest type, areas suited
    to that profile can be highlighted.

====================================================================================

#43  Ordinal Logistic Regression
    file: 43_ordinal_logistic_health_level.xlsx
  >> SCENARIO (narration):
    We predict tree health level (health_level -- ordered: low < medium < high) from age (age), annual precipitation
    (annual_precipitation) and damage score (damage_score). Because the outcome is ordered, ordinal logistic regression
    (cumulative logit) is appropriate: as the damage score rises, how does the probability of shifting to a lower health
    level change, and is precipitation protective?
  >> VARIABLE SELECTION:
    - Dependent (ordered category): health_level
    - Predictors: age, annual_precipitation, damage_score

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Ordinal (Cumulative) Logit   AIC = 407.3   |   ordered: excellent < good < medium < weak
    health_level ~ age + annual_precipitation + damage_score

>> COMMENTARY (narration):
    We modeled tree health level -- an ordered variable -- from age, annual precipitation and damage score. Because
    the outcome is ordered, we used ordinal logistic regression (cumulative logit). This model answers "by how many
    times does the chance of moving to a higher health category change when X rises by one unit". As expected, the
    damage score is the strongest predictor of shifting to a lower health level, while precipitation is protective.
    Under the proportional-odds assumption, a single set of coefficients explains all health thresholds -- which
    simplifies interpretation. In practice, monitoring damage and drought can forecast health deterioration early.

====================================================================================

#44  PLS Regression
    file: 44_pls_spectral_chlorophyll.xlsx
  >> SCENARIO (narration):
    We want to predict leaf chlorophyll content (chlorophyll_ug_g) from 12 spectral bands (wave_1..wave_12). The bands
    are highly correlated (multicollinear) and too numerous relative to the observations; classical regression would be
    unstable. PLS regression solves this by reducing the predictors to components most related to the outcome: how
    accurately does the spectral signature predict chlorophyll?
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: chlorophyll_ug_g
    - Predictors (spectral): wave_1 .. wave_12

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    PLS (2 components)   R²_train = 0.973   R²_CV (5-fold) = 0.968
    chlorophyll_ug_g ~ wave_1 .. wave_12 (spectral)

>> COMMENTARY (narration):
    We wanted to predict leaf chlorophyll content from 12 spectral bands. Because the bands are highly correlated
    and numerous, classical regression would be unstable; PLS regression solved this by reducing the predictors to
    components most related to the outcome. The result is excellent: with just two components, the model explains
    97 percent of the variance in training, and -- crucially -- preserves 96.8 percent in 5-fold cross-validation.
    This shows the model is not memorizing but genuinely generalizable. The practical value is large: leaf spectra
    can estimate chlorophyll -- i.e. photosynthetic capacity and tree health -- non-destructively and quickly. A
    powerful tool for remote sensing and precision forestry.

====================================================================================

#45  Probit Regression (Dose-Response)
    file: 45_probit_herbisit_dose_response.xlsx
  >> SCENARIO (narration):
    We examine the dose-response relationship of a herbicide on a wild weed: against increasing log-dose
    (dose_log_mgL), whether the weed died (died: 0/1). The standard for toxicology and dose-response studies is probit
    regression (normal-CDF link): as the dose rises, how does the probability of death increase, and what is the
    effective dose (e.g. LD50)? Temperature (temperature) can be added as an extra covariate.
  >> VARIABLE SELECTION:
    - Dependent (binary) variable: died
    - Predictor (dose): dose_log_mgL
    - Extra covariate: temperature

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Classification metrics added: accuracy / sensitivity / specificity / AUC (marginal effects + McFadden already shown).

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Probit   AIC = 203.8   Pseudo R² (McFadden) = 0.43   |   marginal effects reported
    died ~ dose_log_mgL + temperature

>> COMMENTARY (narration):
    We examined the dose-response relationship of a herbicide on a wild weed: against increasing log-dose, whether
    the weed died. The standard for toxicology and dose-response studies is probit regression (normal-CDF link). The
    model is strong: McFadden pseudo-R² of 0.43, a very good fit on the logistic/probit scale. Dose, as expected,
    markedly raises the probability of death; temperature also modulates the effect. Probit's most valuable output
    is the effective dose (e.g. LD50 -- the dose affecting half the population): a practical threshold used directly
    to set the correct dosage in herbicide application.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Classification metrics added: accuracy / sensitivity / specificity / AUC (marginal effects + McFadden already shown).

====================================================================================

#46  Tobit Regression (Censored)
    file: 46_tobit_conservation_expenditure_zero_accumulation.xlsx
  >> SCENARIO (narration):
    We model municipalities' forest-conservation expenditure (conservation_expenditure_TL); however, because many
    municipalities spend nothing, the variable piles up at zero (left-censored). Treating the zeros as ordinary numbers
    would bias the estimates. Tobit regression explicitly models the zero pile-up and gives the true effect of budget
    (budget_TL), forest age (forest_age) and awareness score (awareness_score) on spending.
  >> VARIABLE SELECTION:
    - Dependent (left-censored) variable: conservation_expenditure_TL
    - Predictors: budget_TL, forest_age, awareness_score
    - Censoring limit (lower): 0

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Tobit (left-censor = 0)   Log-likelihood = -2631.9   |   coefficient table
    conservation_expenditure_TL ~ budget_TL + forest_age + awareness_score

>> COMMENTARY (narration):
    We modeled municipalities' forest-conservation expenditure; but there was a trap: because many municipalities
    spend nothing, the variable piles up at zero (left-censored). A linear regression treating the zeros as ordinary
    numbers would be biased. Tobit regression solved this by explicitly modeling the zero pile-up. The results show
    the true effect of budget, forest age and awareness on spending: budget and awareness in particular raise
    expenditure. In practice this reveals what drives conservation investment -- such as the multiplier effect of
    investing in awareness on spending, useful for policy design.

====================================================================================

#47  Bayesian Linear Regression
    file: 47_bayesian_lineer_fertilizer_yield.xlsx
  >> SCENARIO (narration):
    We model seedling yield (seedling_yield) from fertilizer doses (nitrogen, phosphorus, potassium) and climate
    (temperature, precipitation); but instead of a classical p-value we want the full probability distribution of each
    coefficient. Bayesian linear regression gives a credible interval and P(B>0) for each effect: e.g. directly
    interpretable statements such as "98% probability that nitrogen increases yield".
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: seedling_yield
    - Predictors: nitrogen_kg_ha, phosphorus_kg_ha, potassium_kg_ha, temperature_C, precipitation_mm

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Bayesian linear (conjugate NIG prior)   per β: posterior mean, 95% credible interval, P(β>0)
    seedling_yield ~ nitrogen + phosphorus + potassium + temperature + precipitation

>> COMMENTARY (narration):
    We modeled seedling yield from fertilizer doses and climate; but instead of a classical p-value we wanted the
    full probability distribution of each effect. Bayesian linear regression gives a credible interval and P(β>0)
    for each coefficient. The beauty is in interpretation: we can make direct, intuitive statements like "98 percent
    probability that nitrogen increases yield" -- without the indirect logic of p-values. The nutrients (especially
    nitrogen) show a strong positive posterior probability on yield. In practice this supports fertilization
    decisions while expressing uncertainty explicitly; the Bayesian approach is especially stable when working with
    little data.

====================================================================================

#48  Nonlinear Regression (Gompertz Growth)
    file: 48_nonlinear_gompertz_age_dbh.xlsx
  >> SCENARIO (narration):
    Tree diameter growth (dbh_cm) with age (age_year) is not linear: fast in youth, then slowing toward an asymptote.
    We cannot capture this S-shaped growth with a linear model; a nonlinear growth curve such as Gompertz is appropriate.
    Nonlinear regression fits the parameters of this form (asymptote, growth rate) to the data: at what age does the
    stand approach maturity?
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: dbh_cm
    - Predictor (time): age_year
    - Model: Gompertz (sigmoid growth)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Gompertz: y = a·exp(-b·exp(-c·x))   R² = 0.9921   RMSE = 1.52
    dbh_cm ~ age_year

>> COMMENTARY (narration):
    Tree diameter growth with age is not linear: fast in youth, then slowing toward an asymptote. We cannot capture
    this S-shaped growth with a linear model, so we fitted a Gompertz curve. The result is extraordinarily good: the
    model explains 99 percent of the variability in diameter (R² = 0.9921), with very low error. The Gompertz
    parameters are biologically meaningful: the asymptote gives the maximum diameter the stand can reach, and the
    growth-rate parameter how quickly it approaches maturity. In practice this curve is very valuable: knowing a
    tree's age we can predict its diameter, or when the stand will approach maturity -- the basis of rotation and
    harvest-timing planning.

====================================================================================

#49  Ridge Regression
    file: 49_ridge_multikolineer_biomass.xlsx
  >> SCENARIO (narration):
    We want to predict biomass (biomass) from 15 soil variables (soil_1..soil_15); but these variables are highly
    correlated (multicollinear), so classical regression coefficients become unstable and inflated. Ridge regression
    (L2 penalty) shrinks large coefficients and removes this instability, giving a more reliable, generalizable model
    under multicollinearity.
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: biomass
    - Predictors: soil_1 .. soil_15
    - Method: Ridge (L2)

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Auto-alpha via cross-validation (RidgeCV): selects the optimal regularization strength automatically.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Ridge (α = 1.0)   R² = 0.508   Adj. R² = 0.468   (n = 200)
    biomass ~ soil_1 .. soil_15 (multicollinear)

>> COMMENTARY (narration):
    We wanted to predict biomass from 15 soil variables; but these are highly correlated (multicollinear), so
    classical regression coefficients would be unstable and inflated. Ridge regression (L2 penalty) removed this
    instability by shrinking large coefficients. R² = 0.51 -- the soil variables explain about half of biomass.
    Ridge's important property: it does not eliminate any variable, keeping them all together to give a balanced,
    generalizable model under multicollinearity. In practice it is the right choice for obtaining a robust estimate
    from many correlated soil measurements; the stability of the model as a whole is foregrounded rather than
    interpreting individual variables.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Auto-alpha via cross-validation (RidgeCV): selects the optimal regularization strength automatically.

====================================================================================

#50  Lasso Regression (Variable Selection)
    file: 50_lasso_30_predictor_5_significant.xlsx
  >> SCENARIO (narration):
    We have 30 candidate predictors (x1..x30) but believe only a few of them truly affect the target. Lasso regression
    (L1 penalty) shrinks the coefficients of unimportant variables exactly to zero; that is, it performs automatic
    variable selection. We thereby reveal the genuinely significant predictors among the 30 in a sparse, interpretable
    model.
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: target
    - Predictors: x1 .. x30
    - Method: Lasso (L1, variable selection)

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Auto-alpha via cross-validation (LassoCV): selects the optimal regularization strength automatically.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Lasso (α = 0.1)   R² = 0.883   Adj. R² = 0.867   |   23 of 30 variables zeroed -> 7 selected
    target ~ x1 .. x30

>> COMMENTARY (narration):
    We had 30 candidate predictors but believed only a few truly affect the target. Lasso regression (L1 penalty)
    did exactly this: it shrank the coefficients of unimportant variables to exactly zero, performing automatic
    variable selection. The result is striking: it eliminated 23 of the 30 variables, leaving 7 genuinely
    significant predictors -- and this sparse model still explained 88 percent of the variance. Lasso's power is
    here: among hundreds of candidates it separates "signal" from "noise" and produces an interpretable, concise
    model. In high-dimensional forestry/environmental data -- dozens of spectral bands, soil parameters -- it is one
    of the most practical ways to find which variables really matter.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Auto-alpha via cross-validation (LassoCV): selects the optimal regularization strength automatically.

====================================================================================

#51  Mediation Analysis
    file: 51_mediation_stress_growth_drying.xlsx
  >> SCENARIO (narration):
    We know that environmental stress (stress) leads to drying in the tree (drying_score); but is the effect direct, or
    does it first slow growth (growth_slowdown) and then cause drying? Mediation analysis splits the total effect of
    stress on drying into direct and indirect (through growth slowdown) components: how much of the stress-drying
    relationship is explained by growth slowdown as a mediator?
  >> VARIABLE SELECTION:
    - Independent (X): stress
    - Mediator (M): growth_slowdown
    - Dependent (Y): drying_score

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Indirect effect (X->M->Y) = 0.3195   SE = 0.035   p < .001   95% CI = (0.255, 0.393)   SIGNIFICANT
    stress -> growth_slowdown -> drying_score

>> COMMENTARY (narration):
    We knew environmental stress leads to drying in the tree; but is the effect direct, or does it first slow growth
    and then cause drying? Mediation analysis solved exactly this. The indirect effect -- stress first slowing
    growth and triggering drying through it -- is 0.32 and statistically significant; the confidence interval lies
    entirely above zero (0.255-0.393). So growth slowdown is a genuine mediator of the stress-drying relationship.
    The practical meaning is important: in combating drying we must target not just stress but stress's effect on
    growth -- e.g. interventions that support growth in stressed trees could break the drying chain early.

====================================================================================

#52  Path / Structural Model
    file: 52_path_structural_health_model.xlsx
  >> SCENARIO (narration):
    We want to test interrelated processes affecting tree health (health_index) in a single model: water status
    (water_status) and nutrient status (nutrient_status) affect the physiological index (physiological_index), which in
    turn determines health. Path analysis estimates these chained/direct relationships simultaneously and shows which
    path is strongest and how well the model fits the data.
  >> VARIABLE SELECTION:
    - Variables: water_status, nutrient_status, physiological_index, health_index
    - Model: path diagram (water/nutrient -> physiological -> health)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    χ² = 24.78 (p < .001)   CFI = 0.926   TLI = 0.882   RMSEA = 0.126
    path: (water, nutrient) -> physiological_index -> health_index

>> COMMENTARY (narration):
    We tested interrelated processes affecting tree health in a single model: water and nutrient status feed the
    physiological index, which in turn determines health. Path analysis estimated these chained relationships at
    once. The fit indices are reasonable: CFI 0.93 reaches the "good" boundary, while RMSEA 0.13 indicates a moderate
    fit. The paths are significant; the proposed chain -- resources first affecting physiology, physiology affecting
    health -- is largely consistent with the data. In practice this shows it is more correct to address tree health
    not directly but through sequential processes: the root of health problems often lies in water/nutrient status
    and physiological capacity.

====================================================================================

#53  Linear Mixed Model (LMM)
    file: 53_lmm_tree_panel_dbh.xlsx
  >> SCENARIO (narration):
    The diameter (dbh_cm) of the same trees (tree_id) was measured over several years (year_t0_t3) and under a treatment
    (treatment). Because the measurements come repeatedly from the same individual, independence is violated and
    classical regression would be biased. The linear mixed model (LMM) takes the tree as a random effect, so it models
    the within-individual correlation and correctly estimates the fixed effect of treatment.
  >> VARIABLE SELECTION:
    - Dependent (continuous) variable: dbh_cm
    - Fixed effects: treatment, year_t0_t3
    - Random effect (grouping): tree_id

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Nakagawa marginal R-squared (fixed effects) and conditional R-squared (fixed + random), beside ICC.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Group = tree_id (60 trees) | Group variance = 7.75   Residual variance = 0.27   ICC = 0.967
    dbh_cm ~ treatment + year_t0_t3 + (1 | tree_id)

>> COMMENTARY (narration):
    The diameter of the same 60 trees was followed over several years, under a treatment. Because the measurements
    come repeatedly from the same individual, independence is violated and classical regression would be biased. The
    LMM took the tree as a random effect and gave a striking result: ICC = 0.97. This means almost all of the
    variability in diameter (97 percent) comes from differences BETWEEN trees, and the within-year measurements of
    the same tree are very consistent. Such a high ICC makes a mixed model mandatory -- ignoring this correlation
    would severely distort the standard errors. By modeling this structure correctly, the model reliably estimated
    the fixed effect of treatment.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Nakagawa marginal R-squared (fixed effects) and conditional R-squared (fixed + random), beside ICC.

====================================================================================

#54  Multiple Imputation (Missing Data)
    file: 54_imputation_missing_inventory.xlsx
  >> SCENARIO (narration):
    In the forest inventory, some trees are missing diameter, height, age or biomass measurements. Deleting incomplete
    rows (listwise) both shrinks the sample and biases the result. Multiple imputation fills the missing values from the
    other variables several times, analyzes each set and pools the results; this gives a complete, unbiased estimate that
    also accounts for uncertainty.
  >> VARIABLE SELECTION:
    - Variables to impute/analyze: dbh_cm, height_m, age, biomass_kg

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    m = 5 imputations | biomass_kg ~ dbh_cm + height_m + age   (missing filled, results pooled)

>> COMMENTARY (narration):
    Some trees in the forest inventory had missing measurements. Deleting incomplete rows (listwise) both shrinks
    the sample and biases the result. Multiple imputation filled the missing values from the other variables five
    times, analyzed each set separately and combined the results. The power of this approach is that it also
    accounts for uncertainty: a single imputation would create false precision, whereas multiple imputation reflects
    the variability among imputations in the standard error. In practice this is the gold standard for producing
    unbiased estimates with correct standard errors from inventories with missing data -- very common in the field.

====================================================================================

#55  GEE (Generalized Estimating Equations)
    file: 55_gee_panel_health_treatment.xlsx
  >> SCENARIO (narration):
    The health (health -- binary: healthy/diseased) of the same trees was followed at successive visits (visit) and under
    a treatment (treatment). With repeated binary outcomes there is within-individual correlation; GEE estimates the
    population-level (marginal) effects, correcting that correlation with a cluster-robust error: how does the treatment
    change the probability of remaining healthy, on average?
  >> VARIABLE SELECTION:
    - Dependent (binary) variable: health
    - Fixed effects: treatment, visit
    - Grouping (subject): tree_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    treatment[T.B]: coef = 0.83, p = 0.050 | health (0/1) ~ treatment + visit, group = tree_id
    DECISION: treatment effect at the significance border

>> COMMENTARY (narration):
    The health (healthy/diseased) of the same trees was followed at successive visits, under a treatment. With
    repeated binary outcomes there is within-individual correlation; GEE estimates the population-level (marginal)
    effects, correcting that correlation with a cluster-robust error. The treatment effect is at the significance
    border (p = 0.050): treatment B raises the chance of staying healthy relative to reference, but the evidence is
    at the threshold. The reason for choosing GEE matters: had we treated the repeated measurements as independent,
    the standard errors would be misleading, creating false precision. The marginal interpretation answers correctly
    "how does treatment change the health probability for an average tree".

====================================================================================

#56  GLMM (Generalized Linear Mixed Model)
    file: 56_glmm_site_panel_insect.xlsx
  >> SCENARIO (narration):
    The insect count (insect_count -- count data) was followed over time (time) at different sites (site_id) as a function
    of precipitation (precipitation). There is both a count distribution and within-site repeated measurement; GLMM uses a
    Poisson link for the count while taking the site as a random effect. It thus correctly estimates the fixed effect of
    precipitation on insect density while modeling site differences.
  >> VARIABLE SELECTION:
    - Dependent (count) variable: insect_count
    - Fixed effects: precipitation, time
    - Random effect (grouping): site_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Poisson GLMM   insect_count ~ precipitation + time + (1 | site_id)
    (count distribution + site random effect modeled)

>> COMMENTARY (narration):
    Insect count was followed over time at different sites, as a function of precipitation. Two challenges combined
    here: a count distribution (non-negative integers) and within-site repeated measurement. GLMM solved both: a
    Poisson link for the count, with the site as a random effect. We thus correctly estimated the fixed effect of
    precipitation on insect density while modeling between-site differences. In practice this lets us answer "how
    does insect count change as precipitation rises" while accounting for each site's own baseline -- the correct
    way to analyze multi-site monitoring studies.

====================================================================================

#57  Regularized Regression (Elastic Net)
    file: 57_regularized_elastic_net.xlsx
  >> SCENARIO (narration):
    We predict the target from 25 features (f1..f25); the variables are both numerous and correlated with one another.
    Lasso alone picks one of a set of correlated variables and drops the others; Ridge eliminates none. Elastic Net
    combines the L1 and L2 penalties, so it both selects variables and keeps correlated groups together; in high-
    dimensional, multicollinear data it is the most balanced choice.
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: target
    - Predictors: f1 .. f25
    - Method: Elastic Net (L1+L2)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Elastic Net   α = 0.080   R²_train = 0.897   R²_CV = 0.861   (target ~ f1..f25)

>> COMMENTARY (narration):
    We predicted the target from 25 features; the variables were both numerous and correlated. Lasso alone picks one
    of a set of correlated variables and drops the others; Ridge eliminates none. Elastic Net combined the L1 and L2
    penalties, so it both selected variables and kept correlated groups together. The result is robust: R² 0.90 in
    training and 0.86 in cross-validation -- the model does not memorize, it generalizes. In high-dimensional,
    multicollinear data Elastic Net is the most balanced choice: it combines Lasso's sparsity with Ridge's stability.
    In environmental and remote-sensing data, working with correlated feature groups, it is an ideal method.

====================================================================================

#58  Robust Regression
    file: 58_robust_regression_outlier_height.xlsx
  >> SCENARIO (narration):
    We model tree height (height_m) by age (age_year); but the data contain a few extreme outliers (e.g. measurement
    error or an abnormal individual). Classical OLS regression is strongly influenced by outliers and the slope becomes
    misleading. Robust regression down-weights the outliers and gives a resilient model faithful to the main trend of the
    data.
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: height_m
    - Predictor: age_year
    - Method: robust (M-estimation)

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Huber-T M-estimator   N = 100 | 8 observations down-weighted as outliers (w < 0.5)
    height_m ~ age_year

>> COMMENTARY (narration):
    We modeled tree height by age; but the data contained a few extreme outliers -- measurement error or an abnormal
    individual. Classical OLS regression is strongly influenced by outliers and the slope becomes misleading. Robust
    regression (Huber-T) down-weighted the outliers and stayed faithful to the main trend of the data: of 100
    observations, 8 were detected as outliers and down-weighted. The practical value is large: measurement errors are
    inevitable in field data; robust methods keep a few erroneous records from distorting the whole model. A large
    difference between the robust and OLS coefficients also reveals how much the outliers affected the classical model
    -- in such cases the robust estimate is more reliable.

====================================================================================

#59  Quantile Regression
    file: 59_quantile_forest_worker_income.xlsx
  >> SCENARIO (narration):
    We model forest workers' income (monthly_income_TL) by education (education_year) and experience (experience_year);
    but we are interested in the tails of the distribution rather than the mean. Classical regression gives only the mean.
    Quantile regression shows how the effect of education changes for low-income (e.g. the 0.10 quantile) and high-income
    (0.90 quantile) workers: is experience more decisive among the low earners?
  >> VARIABLE SELECTION:
    - Dependent (outcome) variable: monthly_income_TL
    - Predictors: education_year, experience_year
    - Quantiles: 0.10, 0.50, 0.90

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Quantile regression: q = 0.10, 0.50, 0.90 | education coefficient varies by quantile (e.g. q=0.10 ~ 0.77)
    monthly_income_TL ~ education_year + experience_year

>> COMMENTARY (narration):
    We modeled forest workers' income by education and experience; but we were interested in the tails of the
    distribution rather than the mean. Classical regression gives only the mean; quantile regression shows how the
    effect changes at different income levels. The results are instructive: the effect of education is not the same
    for low-income (q=0.10) and high-income (q=0.90) workers -- the coefficient varies by quantile. This challenges
    the assumption that "education raises everyone's income equally"; the effect depends on where you are in the
    income distribution. For policy this is very valuable: seeing whether an intervention helps low or high earners
    more is critical for equity-focused decisions.

====================================================================================

#60  ROC Analysis
    file: 60_roc_biomarker_disease.xlsx
  >> SCENARIO (narration):
    We assess how well a biomarker score (biomarker_score) discriminates tree disease (disease: 0/1). ROC analysis plots
    the sensitivity-specificity trade-off at different thresholds; the area under the curve (AUC) summarizes the
    biomarker's discriminative power in a single number (0.5 chance, 1.0 perfect). What is the optimal cut-off, and can
    this biomarker be used for screening?
  >> VARIABLE SELECTION:
    - Score (continuous): biomarker_score
    - Outcome (binary): disease

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    AUC = 0.886 (good discrimination)   Optimal threshold (Youden J) = 0.318   Sensitivity = 0.875
    biomarker_score -> disease (n_+ = 112, n_- = 188)

>> COMMENTARY (narration):
    We assessed how well a biomarker score discriminates tree disease. ROC analysis plotted the sensitivity-
    specificity trade-off at different thresholds; the area under the curve (AUC) summarizes that discriminative
    power in a single number. The result is good: AUC = 0.89 -- the biomarker correctly ranks diseased versus healthy
    trees with 89 percent probability (0.5 chance, 1.0 perfect). The optimal threshold found by Youden J is 0.32,
    with a sensitivity of 87.5 percent there. In practice this shows the biomarker could be used as a screening tool:
    trees above this threshold can be flagged as at-risk and sent for closer examination. The threshold can be tuned
    to the cost of false positives versus false negatives.

====================================================================================

#61  Species Distribution Model (TSS)
    file: 61_tss_type_distribution_model.xlsx
  >> SCENARIO (narration):
    We model the presence/absence of a tree species (type_varligi: 0/1) with environmental variables -- elevation
    (elevation_m), precipitation (precipitation_mm), temperature (temperature_C). The species distribution model predicts
    in which climatic conditions the species is most likely to occur and measures model performance with the TSS (True
    Skill Statistic): in what elevation-climate niche is the species optimal?
  >> VARIABLE SELECTION:
    - Dependent (presence/absence): type_varligi
    - Environmental predictors: elevation_m, precipitation_mm, temperature_C

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    TSS = 0.21 (acceptable)   Sensitivity = 0.62   Specificity = 0.58   N = 500
    type_varligi ~ elevation_m + precipitation_mm + temperature_C

>> COMMENTARY (narration):
    We modeled the presence/absence of a tree species from environmental variables -- elevation, precipitation,
    temperature. The species distribution model predicts in which climatic conditions the species is most likely to
    occur and measures success with the TSS. TSS = 0.21 is "acceptable" but not strong; sensitivity and specificity
    are around 60 percent. This means these three climate variables partly explain the species' distribution but are
    not sufficient alone -- soil, aspect, competition may also be at play. In practice this suggests enriching the
    model with more variables or different algorithms would help; even so, the current form can be used to predict
    the species' coarse climatic niche.

====================================================================================

#62  Confusion Matrix (Classification Evaluation)
    file: 62_confusion_matris_3sinif_diagnosis.xlsx
  >> SCENARIO (narration):
    We evaluate the performance of a diagnostic model across three health classes (actual_label vs prediction_label). The
    confusion matrix cross-tabulates true and predicted classes; from it we compute metrics such as accuracy, precision,
    sensitivity and F1: which classes does the model confuse, and in which class is it weak?
  >> VARIABLE SELECTION:
    - Actual label: actual_label
    - Predicted label: prediction_label

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Accuracy = 0.94   Precision = 0.87   F1 = 0.93   (3 health classes)
    actual_label vs prediction_label

>> COMMENTARY (narration):
    We evaluated the performance of a diagnostic model across three health classes. The confusion matrix
    cross-tabulates the true and predicted classes, from which all metrics are computed. The result is very good:
    accuracy 94 percent, F1 score 0.93. So the model predicts the large majority of classes correctly, with high and
    balanced precision and recall. The real value of the confusion matrix is showing WHICH classes the model
    confuses: if errors concentrate in a particular class despite high overall accuracy, that is where to focus.
    These high scores suggest the model could be used confidently in field diagnosis.

====================================================================================

#63  Random Forest (Classification)
    file: 63_random_forest_forest_type.xlsx
  >> SCENARIO (narration):
    We want to predict a plot's forest type (forest_type) from environmental features -- age, elevation, precipitation,
    soil pH, moisture. Random forest combines hundreds of decision trees to capture nonlinear relationships and
    interactions; it also ranks feature importance: which environmental variable most determines forest type?
  >> VARIABLE SELECTION:
    - Dependent (category): forest_type
    - Predictors: age, elevation, precipitation, pH, moisture

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Test accuracy = 0.838   |   forest_type ~ age + elevation + precipitation + pH + moisture
    (feature importance ranking reported)

>> COMMENTARY (narration):
    We wanted to predict a plot's forest type from environmental features. Random forest combines hundreds of
    decision trees to capture nonlinear relationships and interactions. The test accuracy is 83.8 percent -- with the
    five environmental variables it largely classifies forest type correctly. One of random forest's most valuable
    outputs is feature importance: it shows which environmental variable most determines forest type (usually
    elevation and precipitation stand out). In practice this is used both for mapping (predicting forest type from
    environmental layers) and ecological understanding (which factors are decisive).

====================================================================================

#64  Support Vector Machine (SVM)
    file: 64_svm_spectral_thermal_disease.xlsx
  >> SCENARIO (narration):
    We want to discriminate disease status (disease: 0/1) from spectral (spectral_1, spectral_2) and thermal (thermal_1,
    thermal_2) bands. The classes may not be linearly separable; SVM, via a kernel transformation, forms the widest
    separating margin in a high-dimensional space to capture complex boundaries: how well does the spectral-thermal
    signature separate diseased trees from healthy ones?
  >> VARIABLE SELECTION:
    - Dependent (binary) variable: disease
    - Predictors: spectral_1, spectral_2, thermal_1, thermal_2

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Test accuracy = 0.686   F1 = 0.607   (RBF kernel)   disease ~ spectral_1/2 + thermal_1/2

>> COMMENTARY (narration):
    We tried to discriminate disease status from spectral and thermal bands. The classes may not be linearly
    separable; SVM, via the RBF kernel transformation, forms complex boundaries in a high-dimensional space. The test
    accuracy is 68.6 percent, F1 0.61 -- moderate discrimination. This means the four bands carry part of the disease
    signal but are not sufficient for perfect separation; more bands, better feature engineering, or a different
    kernel could be tried. Still, it is clear the spectral-thermal signature has a connection to disease. Early
    disease detection by remote sensing is promising but requires more work on this dataset.

====================================================================================

#65  Gradient Boosting
    file: 65_gradient_boost_fire_risk.xlsx
  >> SCENARIO (narration):
    We predict fire risk (fire: 0/1) from temperature, moisture, wind (wind), drought index (drought_index) and plant
    density (plant_density). Gradient boosting adds weak trees sequentially, each correcting the previous one's errors, to
    build a high-accuracy model and rank the riskiest factors: which driver most triggers fire?
  >> VARIABLE SELECTION:
    - Dependent (binary) variable: fire
    - Predictors: temperature, moisture, wind, drought_index, plant_density

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Test accuracy = 0.55   F1 = 0.58   |   fire ~ temperature + moisture + wind + drought_index + plant_density
    (feature importance ranking reported)

>> COMMENTARY (narration):
    We tried to predict fire risk from five environmental drivers. Gradient boosting adds weak trees sequentially,
    each correcting the previous one's errors, to build strong models. But on this dataset the test accuracy is 55
    percent -- slightly above chance, but low. This is an honest result: fire, by nature, depends on many interacting
    and stochastic factors, hard to predict from five variables. The model's feature importance is still valuable,
    highlighting the riskiest drivers (usually drought and temperature). The practical message: fire-risk modeling
    needs richer data (moisture history, topography, fuel load) and more observations; the current model is a
    starting point.

====================================================================================

#66  K-Means Clustering
    file: 66_kmeans_forest_type_morfometri.xlsx
  >> SCENARIO (narration):
    Without using labels, we want to split trees into natural groups by their morphometric similarity -- diameter (dbh_cm),
    height (height_m), biomass (biomass_kg). K-means assigns trees to the nearest cluster centers, bringing similar
    individuals together: how many distinct "growth types" are there, and which morphological profiles separate the
    clusters?
  >> VARIABLE SELECTION:
    - Clustering variables: dbh_cm, height_m, biomass_kg

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    3 clusters | mean dbh/height/biomass reported per cluster | inertia (elbow) curve
    cols: dbh_cm, height_m, biomass_kg

>> COMMENTARY (narration):
    Without using labels, we split trees into three natural clusters by their morphometric similarity -- diameter,
    height, biomass. K-means assigned similar individuals to the nearest cluster center. The cluster summaries show
    each cluster's profile: typically one cluster gathers small/young trees, one medium, one large/mature -- so
    "growth types" or maturity classes emerge. The elbow (inertia) curve confirms three clusters is a reasonable
    choice. In practice this forms the basis for splitting a stand into subgroups from unlabeled data -- e.g.
    differentiating silvicultural interventions by maturity class.

====================================================================================

#67  Hierarchical Clustering
    file: 67_hierarchical_type_protein_profili.xlsx
  >> SCENARIO (narration):
    We want to examine how samples cluster by expression profiles made of ten protein measurements (protein_1..protein_10).
    Hierarchical clustering merges samples step by step by similarity and builds a dendrogram; we cut the tree to read the
    natural groups without specifying the number of clusters in advance: do the protein profiles separate by species (type)?
  >> VARIABLE SELECTION:
    - Clustering variables: protein_1 .. protein_10

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    3 clusters (Ward linkage) | 10 protein profiles (protein_1..protein_10) | dendrogram

>> COMMENTARY (narration):
    We examined how samples cluster by expression profiles of ten protein measurements. Hierarchical clustering
    merged samples step by step by similarity and built a dendrogram; we cut the tree into three clusters. Unlike
    k-means, it shows the similarity hierarchy without imposing the number of clusters in advance -- you can look at
    the dendrogram and choose the natural cut. The protein profiles separating into three groups shows the samples
    differ biochemically, probably corresponding to species or physiological states. In genomic/proteomic forest
    studies this is the standard method for exploratory grouping of samples.

====================================================================================

#68  DBSCAN (Spatial Clustering)
    file: 68_dbscan_fire_hotspot.xlsx
  >> SCENARIO (narration):
    We want to find concentration points (hotspots) of fire events (event_id) by their geographic location (latitude,
    longitude). DBSCAN performs density-based clustering: nearby events form a cluster, isolated events are separated as
    "noise". There is no need to specify the number of clusters in advance; it directly reveals irregularly shaped fire
    hotspots and outlier events.
  >> VARIABLE SELECTION:
    - Coordinates: latitude, longitude

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Density-based clustering (lat/lon) | clusters + noise separated
    cols: latitude, longitude

>> COMMENTARY (narration):
    We wanted to find concentration points of fire events by their geographic location. DBSCAN performs
    density-based clustering: nearby events form a cluster, isolated events are separated as "noise". It differs from
    k-means in two important ways: no need to specify the number of clusters in advance, and it can capture
    irregularly shaped clusters. We thus directly distinguished fire hotspots -- regions where events concentrate --
    and a few isolated events. In practice this is directly usable for directing intervention resources and
    prevention efforts to the densest risk areas; it is one of the most practical tools for spatial hotspot detection.

====================================================================================

#69  PCA (Principal Component Analysis)
    file: 69_pca_forest_measurements_6_feature.xlsx
  >> SCENARIO (narration):
    We want to reduce six correlated forest measurements (feature_1..feature_6) to fewer independent dimensions. Principal
    component analysis (PCA) forms new axes (components) that explain the most shared variance among the variables; the
    first few components usually carry the bulk of the variance. We thus summarize the data, visualize it and remove
    multicollinearity.
  >> VARIABLE SELECTION:
    - Variables: feature_1 .. feature_6

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    KMO = 0.769 (good) | Explained variance: PC1 47.0%, PC2 18.2%, PC3 16.3% (first 3 ~ 81.5%)
    cols: feature_1 .. feature_6

>> COMMENTARY (narration):
    We wanted to reduce six correlated forest measurements to fewer independent dimensions. First the KMO value
    (0.77) confirmed the data are suitable for PCA. The result is efficient: the first component alone carries 47
    percent of the variance, the first three 81.5 percent. So we can summarize most of the information in six
    variables with three dimensions. PC1 usually behaves like a "general size" axis -- gathering the shared variance
    of correlated measurements. The practical benefit: in later analyses we use a few components instead of six
    variables, simplifying the data and removing multicollinearity; we can also visualize samples on a 2D score plot.

====================================================================================

#70  t-SNE (Nonlinear Dimensionality Reduction)
    file: 70_tsne_leaf_features.xlsx
  >> SCENARIO (narration):
    We want to visualize high-dimensional data of twelve leaf features (feature_1..feature_12) in two dimensions; the goal
    is to see whether the species (type) separate naturally. t-SNE produces a nonlinear embedding that preserves local
    neighborhoods: similar leaves appear close and different species as separate clusters. A powerful visual for
    exploratory classification.
  >> VARIABLE SELECTION:
    - Features: feature_1 .. feature_12
    - Color/label (exploration): type

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    KL Divergence = 0.0065 (very good fit) | 12 features -> 2D embedding (color: type)
    cols: feature_1 .. feature_12

>> COMMENTARY (narration):
    We visualized high-dimensional data of twelve leaf features in two dimensions; the goal was to see whether the
    species separate naturally. t-SNE produces a nonlinear embedding that preserves local neighborhoods: similar
    leaves appear close and different species as separate clusters. With a KL divergence of 0.0065 the fit is very
    good -- the embedding faithfully preserved the original high-dimensional structure. An important note: on a t-SNE
    map the axis values are not interpreted; only the relative positions of points are meaningful. In practice this
    is a powerful exploratory visual for seeing whether species separate morphometrically -- the "get to know the
    data" step before classification.

====================================================================================

#71  MDS (Multidimensional Scaling)
    file: 71_mds_leaf_morfometri.xlsx
  >> SCENARIO (narration):
    Based on leaf morphometry -- length, width, area, petiole, pointedness -- we want to show the similarity of samples on a
    two-dimensional map. Multidimensional scaling (MDS) produces a low-dimensional layout that preserves the distances
    between samples: do the natural groups (group_natural) separate on the map, and which leaves are morphologically close?
  >> VARIABLE SELECTION:
    - Measurement variables: length_cm, width_cm, area_cm2, petiole_cm, three_sivrilik
    - Label (exploration): group_natural

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Stress (Kruskal-1) = 0.066 (acceptable/good) | 5 measurements -> 2D map
    length, width, area, petiole, pointedness

>> COMMENTARY (narration):
    Based on leaf morphometry -- length, width, area, petiole, pointedness -- we showed the similarity of samples on
    a two-dimensional map. MDS produces a low-dimensional layout that preserves the distances between samples. The
    quality measure, Kruskal stress, is 0.066 -- very good; the five-dimensional distance structure is represented in
    two dimensions with almost no distortion. In practice this map shows at a glance which leaves are morphologically
    close and whether natural groups separate -- an intuitive tool for taxonomic and phenotypic exploration.

====================================================================================

#72  UMAP (Dimensionality Reduction)
    file: 72_umap_spectral_lower_type.xlsx
  >> SCENARIO (narration):
    We want to embed data of twenty-five spectral channels (channel_1..channel_25) in two dimensions to see whether the
    sub-types (lower_type) separate. UMAP performs a nonlinear reduction that preserves both local and global structure and
    is generally faster than t-SNE: do the spectral signatures arrange the sub-types into separate clusters?
  >> VARIABLE SELECTION:
    - Features (spectral): channel_1 .. channel_25
    - Label (exploration): lower_type

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    UMAP (2D) | 25 spectral channels -> 2D embedding (color: lower_type)

>> COMMENTARY (narration):
    We embedded data of twenty-five spectral channels in two dimensions to see whether the sub-types separate. UMAP
    is a nonlinear reduction that preserves both local and global structure and is generally faster than t-SNE. Its
    advantage over t-SNE is preserving GLOBAL distances between clusters better -- so not only "which points are
    close" but "how far apart clusters are" becomes meaningful. The spectral signatures arranging the sub-types into
    separate clusters is exploratory evidence that species/sub-type discrimination by remote sensing is feasible; a
    powerful way to understand the data's structure before building a classification model.

====================================================================================

#73  Cronbach's Alpha (Reliability)
    file: 73_cronbach_21madde_satisfaction.xlsx
  >> SCENARIO (narration):
    We measure the internal consistency of the 21 items (M01..M21) of a forest-recreation satisfaction scale. Cronbach's
    alpha shows how consistently the items measure the same construct (satisfaction): 0.70+ acceptable, 0.80+ good. The
    item-total correlations and the "alpha if item deleted" table also reveal items that weaken the scale.
  >> VARIABLE SELECTION:
    - Scale items: M01 .. M21

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    α (Cronbach) = 0.942 (excellent)   mean inter-item r = 0.436   (21 items)

>> COMMENTARY (narration):
    We measured the internal consistency of the 21 items of a forest-recreation satisfaction scale. Cronbach's alpha
    shows how consistently the items measure the same construct -- satisfaction. The result is excellent: alpha =
    0.94. On the scale 0.70 is acceptable and 0.80 is good; 0.94 shows the items measure the same concept with very
    high internal consistency. The mean inter-item correlation of 0.44 is healthy -- neither too low (unrelated items)
    nor too high (redundant repetition). In practice this confirms the scale is reliable; the sum/mean of the items
    can be confidently used as a satisfaction score. Very high alpha can also hint at item redundancy; the "alpha if
    item deleted" table shows simplification opportunities.

====================================================================================

#74  Likert Analysis
    file: 74_likert_15madde_3boyut.xlsx
  >> SCENARIO (narration):
    We analyze the responses of a three-dimensional (quality, service, price) 15-item Likert scale. Likert analysis
    summarizes the item distributions, means and dimension scores of each dimension and visualizes the response pattern: in
    which dimension are visitors more satisfied, and which items score low?
  >> VARIABLE SELECTION:
    - Quality dimension: quality_M1 .. quality_M5
    - Service dimension: service_M1 .. service_M5
    - Price dimension: price_M1 .. price_M5

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    5-point Likert (1-5) | 3 dimensions: quality (5), service (5), price (5) | dimension scores + distributions

>> COMMENTARY (narration):
    We analyzed the responses of a three-dimensional -- quality, service, price -- 15-item Likert scale. Likert
    analysis summarizes each dimension's item distributions, means and dimension scores, and visualizes the response
    pattern (how many said "strongly agree", etc.). This way we capture dimensional differences that a single overall
    average would hide: e.g. visitors may be satisfied with quality but find price high. In practice this directly
    shows where to prioritize in service improvement -- the lowest-scoring dimension and items are clear targets for
    intervention.

====================================================================================

#75  Exploratory Factor Analysis (EFA)
    file: 75_efa_18madde_3faktor.xlsx
  >> SCENARIO (narration):
    We do not know how many latent dimensions (factors) lie behind eighteen items (M01..M18) and want to discover them.
    Exploratory factor analysis (EFA) examines the shared variance among items to group them into a few factors; the factor
    loadings show which item measures which dimension. KMO and Bartlett's test also check the suitability of the data.
  >> VARIABLE SELECTION:
    - Scale items: M01 .. M18

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    KMO = 0.897 (very good)   Bartlett χ²(153) = 3298.70, p < .001   |   18 items
    -> factor analysis suitable

>> COMMENTARY (narration):
    We did not know how many latent dimensions lie behind eighteen items and wanted to discover them. First we ran
    two checks: the KMO sampling adequacy of 0.90 (very good) and Bartlett's sphericity test p < .001 -- both say the
    data are highly suitable for factor analysis. Doing EFA without this would be meaningless. Exploratory factor
    analysis examined the shared variance among items to group them into a few factors; the factor loadings show
    which item measures which dimension. In practice this is the way to DISCOVER a scale's structure from the data --
    revealing how many sub-dimensions the items actually split into and which item belongs where.

====================================================================================

#76  ICC (Intraclass Correlation)
    file: 76_icc_3uzman_tree_height.xlsx
  >> SCENARIO (narration):
    We assess the consistency of three experts (expert) independently measuring the height (height_m) of the same trees. The
    intraclass correlation coefficient (ICC) gives the inter-rater reliability: ICC 0.75+ good, 0.90+ excellent agreement.
    Are the measurements independent of the rater, or does the result change depending on who measures?
  >> VARIABLE SELECTION:
    - Measurement variable: height_m
    - Rater (grouping): expert
    - Object: tree_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    ICC(1,1) = 0.980   ICC(A,1) = 0.980 (single rater) | ICC(.,k) ~ 0.993 (average)   p < .001   -> Excellent
    3 experts, height_m

>> COMMENTARY (narration):
    We assessed the consistency of three experts independently measuring the height of the same trees. The
    intraclass correlation coefficient (ICC) gives the inter-rater reliability. The result is nearly perfect: ICC =
    0.98 for a single expert's measurement, 0.99 for the average of experts -- both "excellent". So the measurements
    are independent of the rater; whoever measures, the result is practically the same. This is an important
    assurance: the height-measurement protocol is highly reliable, and data from different teams can be safely
    combined. ICC's six types answer different questions (absolute agreement vs consistency; single vs average); here
    all are excellent, presenting a strong agreement picture.

====================================================================================

#77  Confirmatory Factor Analysis (CFA)
    file: 77_cfa_12madde_3faktor_dogrulayici.xlsx
  >> SCENARIO (narration):
    We test the assumption that twelve items load onto three pre-defined factors -- management (management), worker (worker)
    and site (site). Whereas EFA discovers dimensions, confirmatory factor analysis (CFA) tests whether this theoretically
    specified structure fits the data; the scale structure is confirmed with fit indices (CFI, RMSEA) and factor loadings.
  >> VARIABLE SELECTION:
    - Factor 1 (management): management_M1 .. management_M4
    - Factor 2 (worker): worker_M1 .. worker_M4
    - Factor 3 (site): site_M1 .. site_M4

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    χ² = 52.40, p = 0.42 (ns = good fit)   CFI = 0.9993   |   3 factors: management, worker, site
    DECISION: model fits the data well

>> COMMENTARY (narration):
    We tested the assumption that twelve items load onto three pre-defined factors -- management, worker, site.
    Whereas EFA DISCOVERS dimensions, CFA TESTS whether this theoretically specified structure fits the data. The
    result is very good: chi-square is non-significant (p = 0.42) -- which in CFA is what you want, meaning no
    significant departure between model and data. With CFI = 0.999 the fit is nearly perfect. So the three-factor
    structure is confirmed: the 12 items genuinely separate cleanly into management, worker and site dimensions. In
    practice this proves the scale's construct validity -- the scale really measures the three concepts it claims to.

====================================================================================

#78  Survey Means (Weighted Mean)
    file: 78_survey_means_forest_worker_geliri.xlsx
  >> SCENARIO (narration):
    We estimate the mean income (monthly_income) of forest workers; but the sample is not simple random -- it has a stratified
    (by region), weighted (weight), clustered (psu_province) design. The complex-sample mean accounts for these design weights
    and clustering to give the population mean and the correct standard error; ignoring the weights would lead to a biased
    estimate.
  >> VARIABLE SELECTION:
    - Outcome (continuous): monthly_income
    - Weight: weight
    - Stratum: region
    - Cluster (PSU): psu_province

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Weighted mean income M̂ = 4,876 TL   SE = 33.5   95% CI = [4,809, 4,943]   CV = 0.69%   (Taylor SE)
    stratum = region, cluster = psu_province

>> COMMENTARY (narration):
    We estimated the mean income of forest workers; but the sample was not simple random -- it had a stratified
    (by region), weighted, clustered design. The complex-sample mean accounted for these design weights and
    clustering to give the population mean and the correct standard error (Taylor linearization): a mean income of
    4,876 TL with a very narrow confidence interval. Ignoring the weights would have biased the estimate. In practice
    this is the standard for producing population-GENERALIZABLE estimates with correct uncertainty from survey data
    -- indispensable in official statistics and field surveys.

====================================================================================

#79  Survey Frequency (Weighted Proportion)
    file: 79_survey_freq_insured_ratio.xlsx
  >> SCENARIO (narration):
    We estimate the rate of being insured (insured: 0/1) among forest workers under a complex sample. Weighted frequency
    analysis uses the design weights (weight), strata and clusters (psu) to correctly compute the population proportion and its
    confidence interval: what is the insurance rate, and how does it vary across regions?
  >> VARIABLE SELECTION:
    - Outcome (categorical/binary): insured
    - Weight: weight
    - Stratum: region
    - Cluster (PSU): psu

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Insured proportion p̂ = 0.645 (64.5%)   SE = 0.036   95% CI = [0.570, 0.714]   (Taylor SE)
    stratum = region, cluster = psu

>> COMMENTARY (narration):
    We estimated the rate of being insured among forest workers under a complex sample. Weighted frequency analysis
    used the design weights, strata and clusters to compute the population proportion and its correct confidence
    interval: about 64.5 percent of workers are insured, with a 57-71 percent interval. The critical point is that
    the standard error is computed per the design: because of the cluster effect, the true uncertainty is larger than
    a simple-random-sample assumption would give. In practice this proportion and interval are directly usable for
    social policy (intervention toward uninsured workers) -- a population-generalizable estimate.

====================================================================================

#80  Survey Total (Weighted Total)
    file: 80_survey_total_forest_area_tabakali.xlsx
  >> SCENARIO (narration):
    From a stratified forest inventory we want to estimate the total cultivated area (cultivated_area_ha) and total production.
    Weighted total analysis uses each forest's design weight (weight) and stratum (stratum) to give the total estimate that
    generalizes from the sample to the whole population, along with its standard error: what is the total forest area/production?
  >> VARIABLE SELECTION:
    - To total (continuous): cultivated_area_ha
    - Weight: weight
    - Stratum: stratum

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Total cultivated area T̂ = 969,817 ha   SE = 16,367   95% CI = [937,655, 1,001,979]   (Taylor SE)
    stratum = stratum

>> COMMENTARY (narration):
    From a stratified forest inventory we estimated the total cultivated area. Weighted total analysis used each
    forest's design weight and stratum to give the total estimate generalizing from the sample to the whole
    population, with its standard error: about 970 thousand hectares, with a 938 thousand - 1.002 million hectare
    interval. The TOTAL estimate, not the mean, is critical in resource inventory and planning -- it answers "how
    much total forest/production is there". The stratified design improves the precision of the estimate; the
    standard error clearly states how reliable the total is. Carbon stock, harvest quota and land-management decisions
    rest directly on such total estimates.

====================================================================================

#81  Survey Regression (Weighted Linear)
    file: 81_survey_reg_forest_worker_geliri.xlsx
  >> SCENARIO (narration):
    We model forest workers' income (monthly_income) by age, education (education_year) and sex; but the data have a complex
    sample design (weight, stratum=region, cluster=psu). Weighted linear regression accounts for the design weights and
    clustering and produces cluster-robust standard errors: what is the effect of education on income, generalizable to the
    population?
  >> VARIABLE SELECTION:
    - Dependent (outcome): monthly_income
    - Predictors: age, education_year, sex
    - Weight: weight | Stratum: region | Cluster: psu

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    R² = 0.297   education_year: coef = 271.5, p < .001 (95% CI [205, 338])   (stratum=region, cluster=psu)
    monthly_income ~ age + education_year

>> COMMENTARY (narration):
    We modeled forest workers' income by age and education; the data had a complex sample design (weight, stratum,
    cluster). Weighted linear regression accounted for the design weights and clustering, producing cluster-robust
    standard errors. The result is clear: each additional year of education raises monthly income by about 271 TL on
    average, and this effect is highly significant (p < .001), with a wholly positive confidence interval. The model
    explains about 30 percent of the variability in income. Importantly, the estimate is population-generalizable:
    thanks to the design weights it represents the entire worker population, not the sample. This strong, unbiased
    effect of education on income is a concrete basis for vocational-education policy.

====================================================================================

#82  Survey Logistic (Weighted Logistic)
    file: 82_survey_logistic_is_accident.xlsx
  >> SCENARIO (narration):
    We model the probability of having a work accident (is_accident: 0/1) under a complex sample; the predictors are age,
    experience and being insured (insured). Weighted logistic regression uses the design weights and clustering (psu) to give
    odds ratios and cluster-robust confidence intervals: does accident probability decrease significantly as experience rises?
  >> VARIABLE SELECTION:
    - Dependent (binary): is_accident
    - Predictors: age, experience, insured
    - Weight: weight | Stratum: region | Cluster: psu

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    experience: coef = -0.070, OR = 0.932, p < .001 (95% CI [0.90, 0.97])   (stratum=region, cluster=psu)
    is_accident ~ age + experience + insured

>> COMMENTARY (narration):
    We modeled the probability of having a work accident under a complex sample. The most striking finding is in
    experience: each additional year lowers the accident odds to 0.93 times (OR = 0.932, p < .001) -- so as
    experience rises, accident risk falls significantly and consistently. The confidence interval lies entirely below
    1, so the protective effect is robust. Because weighted logistic regression gives these odds ratios with design
    weights and cluster-robust errors, the result is population-generalizable. A very valuable occupational-safety
    finding in practice: inexperienced workers' risk is markedly higher, so training and supervision for newcomers
    yield the highest return in reducing accidents.

====================================================================================

#83  GAM (Generalized Additive Model)
    file: 83_gam_pm25_temperature_u_shape.xlsx
  >> SCENARIO (narration):
    We relate air quality (PM25) to temperature (temperature_C); but the relationship is not linear -- PM25 may rise at both
    very low and very high temperatures in a U-shape. GAM applies smooth functions to the predictors to capture such nonlinear
    patterns without specifying the form in advance: what is the true shape of the temperature-PM25 relationship, and where is
    the threshold?
  >> VARIABLE SELECTION:
    - Dependent (outcome): PM25
    - Smooth predictors: temperature_C, moisture_pct

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Pseudo R² = 0.75   McFadden adj = 0.30   edf = 11.28   AIC = 2364   (N = 365)
    PM25 ~ smooth(temperature_C) + smooth(moisture_pct)

>> COMMENTARY (narration):
    We related air quality to temperature; but the relationship is not linear -- PM25 may rise at both very low and
    very high temperatures in a U-shape. GAM captured such nonlinear patterns by applying smooth functions to the
    predictors, without specifying the form in advance. The model is strong: 75 percent explained variance. The
    effective degrees of freedom (edf) of 11.3 indicate the relationship departs markedly from linear, carrying a
    curved structure (edf = 1 would be linear). In practice GAM's value is here: instead of a simple "temperature
    raises PM25", it reveals the true shape of the relationship -- where it dips, where it peaks. Ideal for capturing
    nonlinear thresholds in environmental and climate data.

====================================================================================

#84  Discriminant Analysis (LDA)
    file: 84_diskriminant_3tur_5ozellik.xlsx
  >> SCENARIO (narration):
    We want to find the linear combinations that discriminate three tree species (type) from five morphological features --
    diameter (dbh), height (height), biomass, crown, bark. Discriminant analysis (LDA) forms the axes that best separate the
    classes and gives the classification accuracy and cross-validated performance: how accurately can the species be separated
    with these five measurements, and which feature is most discriminating?
  >> VARIABLE SELECTION:
    - Dependent (category): type
    - Predictors: dbh, height, biomass, crown, bark

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Accuracy = 0.911   CV accuracy (5-fold) = 0.917 ± 0.030   LD1 carries 98.8% of variance
    type ~ dbh + height + biomass + crown + bark

>> COMMENTARY (narration):
    We found the linear combinations that discriminate three tree species from five morphological features.
    Discriminant analysis (LDA) forms the axes that best separate the classes. The result is very good: classification
    accuracy 91 percent, and 92 percent in cross-validation -- the model does not memorize, it generalizes. The most
    striking detail is that the first discriminant axis (LD1) carries 98.8 percent of the variance: the species
    actually separate along a single dominant axis. This shows the five measurements separate species largely via a
    common "size" dimension. In practice this proves we can classify species with high accuracy from simple
    morphological measurements -- e.g. speeding up field identification.

====================================================================================

#85  Conditional Logit (Choice Model)
    file: 85_conditional_logit_transport_choice.xlsx
  >> SCENARIO (narration):
    To carry forest products to market, each carrier (person_id) chooses one of several transport options (mode); for each
    chooser one option is selected (selected: 0/1). Each option's duration (duration_minute), fee (fee_TL) and distance
    (distance) differ. Conditional logit does what classical logistic cannot: it models a choice made among several options
    based on the option attributes: as fee and duration rise, how does the probability of choosing that option fall?
  >> VARIABLE SELECTION:
    - Chooser: person_id
    - Chosen (0/1): selected
    - Explanatory: duration_minute, fee_TL, distance

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    McFadden R² = 0.21   |   duration: β = -0.029 (p = .13)   fee: β = -0.004 (p = .41)
    chooser = person_id, chosen = selected, n_chooser = 200 (4 alternatives/chooser)

>> COMMENTARY (narration):
    We modeled each carrier choosing one of several transport options to take forest products to market. Conditional
    logit does what classical logistic cannot: it models a choice made among several options based on the option
    attributes (duration, fee, distance). The model shows a reasonable fit overall: McFadden R² = 0.21, which on this
    scale is already "good fit" (0.2-0.4; not to be confused with classical R²). Both coefficients are negative -- as
    duration and fee rise, the probability of choosing that option falls, which is intuitively correct. However, in
    this sample the coefficients do not individually reach statistical significance; the direction is as expected, but
    a larger sample may be needed for firm evidence. Even so, the model shows the choice behavior follows a rational
    pattern.

====================================================================================

#86  Kaplan-Meier (Survival Curve)
    file: 86_kaplan_meier_street_agaci.xlsx
  >> SCENARIO (narration):
    We follow how long street trees survive after planting: follow-up time (followup_year) and death status (tree_died: 0/1).
    Some trees are still alive at the end of follow-up (censored). Kaplan-Meier handles this censored data correctly to plot the
    survival curve against time; curves for groups such as irrigation (irrigation) are compared with the log-rank test: do
    irrigated trees live longer?
  >> VARIABLE SELECTION:
    - Time: followup_year
    - Event (1=died, 0=censored): tree_died
    - Grouping: irrigation

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Median survival: not reached (>50% alive at end of follow-up) | Log-rank χ²(1) = 25.30, p < .001
    time = followup_year, event = tree_died, group = irrigation (n = 200)

>> COMMENTARY (narration):
    We followed how long street trees survive after planting and compared irrigated versus non-irrigated groups.
    Kaplan-Meier handled the censored data (trees still alive at the end of follow-up) correctly to plot the survival
    curve. Two striking results: first, the median survival was not reached -- so more than half the trees are still
    alive during follow-up, a good survival sign. Second and most important, the log-rank test: chi-square 25.30, p
    below one in a thousand. So the survival curves of the irrigation groups differ significantly; irrigated trees
    live markedly longer. The practical message is clear: irrigation strongly increases survival in street trees -- it
    should be an investment priority in urban afforestation.

====================================================================================

#87  Cox Regression
    file: 87_cox_regression_tree_death.xlsx
  >> SCENARIO (narration):
    We look for the continuous factors affecting tree death risk (tree_died): age, stress score (stress_score), temperature,
    drought (drought). The Cox proportional-hazards model gives, in censored survival data, each predictor's effect on the death
    hazard as a hazard ratio (HR): by how many times does death risk rise as the stress score increases, and which factor is most
    decisive?
  >> VARIABLE SELECTION:
    - Time: followup_year
    - Event: tree_died
    - Predictors: age, stress_score, temperature, drought

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Proportional-hazards assumption test (Schoenfeld): per-covariate chi-square/p; significant = PH violated, consider a time-varying effect.

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    stress_score: HR = 1.215, p < .001 (95% CI [1.15, 1.28]) | C-index = 0.689
    time = followup_year, event = tree_died ~ age + stress_score + temperature + drought

>> COMMENTARY (narration):
    We used the Cox proportional-hazards model to find the continuous factors affecting tree death risk. The most
    prominent factor is the stress score: hazard ratio 1.215, p below one in a thousand. This means that when the
    stress score rises by one unit, the death hazard rises by 21.5 percent -- the confidence interval lies entirely
    above 1, so it is a robust risk factor. The model's discriminative power, C-index 0.69, is moderate: stress,
    temperature, drought and age partly predict death risk. Cox's strength is giving each factor's effect as a hazard
    ratio in censored survival data. In practice this shows stressed trees should be prioritized for monitoring, and
    interventions that reduce stress could directly lower death risk.

  >> ADVANCED PARAMETERS (optional in the form — what they do):
    - Proportional-hazards assumption test (Schoenfeld): per-covariate chi-square/p; significant = PH violated, consider a time-varying effect.

====================================================================================

#88  AFT (Accelerated Failure Time / Weibull)
    file: 88_aft_weibull_tree_omru.xlsx
  >> SCENARIO (narration):
    We want to model tree lifespan directly: do rearing (rearing) and species (type) lengthen or shorten the survival time?
    Whereas Cox focuses on the hazard ratio, the AFT (Weibull) model gives how much the event is "accelerated/decelerated" in
    time units: by what percentage does a given rearing method extend tree lifespan? The parametric survival interpretation is
    more intuitive.
  >> VARIABLE SELECTION:
    - Time: age_year
    - Event (1=died): dead
    - Predictors: rearing, type

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Weibull AFT | type_B: coef = -0.432, time ratio = 0.649, p < .001
    time = age_year, event = dead ~ type

>> COMMENTARY (narration):
    We modeled tree lifespan directly: does species lengthen or shorten survival time? Whereas Cox focuses on the
    hazard ratio, the AFT (Weibull) model gives how much the event is accelerated/decelerated in time units -- a more
    intuitive interpretation. The result is significant: species B's time ratio is 0.649 (p < .001). So B trees live
    about 35 percent shorter than the reference species; death is "accelerated" in them. The confidence interval lies
    entirely below 1, a robust effect. In practice this shows species choice determines not just growth but LIFESPAN;
    long-lived species should be preferred, especially in long-rotation or permanent landscape plantings. Parametric
    survival is very useful in practice because it gives directly interpretable statements like "lives x percent
    longer/shorter".

====================================================================================

#89  Competing Risks
    file: 89_competing_risks_3olum_cause.xlsx
  >> SCENARIO (narration):
    Trees can die from more than one cause (event_status: 0=censored, 1/2/3=different causes of death). Standard Cox treats one
    cause as censored and gives biased estimates. Competing-risks analysis models each cause of death -- with the others as
    competitors -- separately (cause-specific Cox + Cumulative Incidence): which cause of death does age/stress particularly
    increase?
  >> VARIABLE SELECTION:
    - Time: followup_year
    - Event type (0=censored, 1/2/3=cause): event_status
    - Predictors: age, stress_score

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Event types: 1 (n=124), 2 (n=209), 3 (n=7) | censored = 60 | cause-specific Cox + CIF
    time = followup_year, event type = event_status ~ age + stress_score

>> COMMENTARY (narration):
    Trees can die from more than one cause; this data has three causes of death and censored (still alive)
    observations. Standard Cox would treat one cause as censored and give biased estimates. Competing-risks analysis
    models each cause of death -- with the others as competitors -- separately (cause-specific Cox) and computes
    cumulative incidence functions. This lets us answer "which cause of death does age and stress particularly
    increase" -- far more informative than "does it raise death risk". The most frequent cause is the second type
    (209 events); the third is rare (7). In an academic report both cause-specific hazard ratios and CIF should be
    presented; they answer different questions -- one "which is the risk factor", the other "what is the actual
    probability".

====================================================================================

#90  Time-Dependent Cox
    file: 90_td_cox_stress_zamanli.xlsx
  >> SCENARIO (narration):
    A tree's stress level (stress) changes over time; using a single baseline value would be misleading. Time-dependent Cox
    structures the data in long-format with (t_start, t_stop] intervals and uses the current stress value in each interval: how
    does the stress at that moment affect the death risk at that moment? Time-varying interventions such as a drug (drug) are thus
    modeled correctly too.
  >> VARIABLE SELECTION:
    - Start/stop: t_start, t_stop
    - Event: event
    - Time-dependent predictors: stress, drug
    - Subject: tree_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    stress: HR = 1.535, p < .001 (95% CI [1.26, 1.88]) | log-likelihood = -244.4
    (t_start, t_stop] long-format, event = event ~ stress + drug, subject = tree_id

>> COMMENTARY (narration):
    A tree's stress level changes over time; using a single baseline value would be misleading. Time-dependent Cox
    structures the data in (t_start, t_stop] intervals and uses the CURRENT stress value in each interval. The result
    is striking: the instantaneous stress strongly raises the death hazard -- HR = 1.535, p below one in a thousand.
    So when the stress at that moment rises by one unit, the death hazard rises by 53.5 percent; this is a stronger,
    more realistic effect than a classical (fixed) Cox could capture, because it accounts for the fluctuation of
    stress over time. In practice this is very important: in tree-health monitoring, the CURRENT stress state -- not a
    past measurement -- is decisive for risk. Time-varying treatments/interventions are thus modeled correctly too.

====================================================================================

#91  Survey-PHREG (Weighted Cox)
    file: 91_survey_phreg_weighted_cox.xlsx
  >> SCENARIO (narration):
    We run the survival analysis under a complex sample: the trees come from a weighted (weight) and clustered (psu) design.
    Survey-PHREG builds a weighted Cox regression and gives cluster-robust standard errors; ignoring the design would lead to
    biased hazard ratios: what is the effect of the stress score on death risk, generalizable to the population?
  >> VARIABLE SELECTION:
    - Time: followup
    - Event: died
    - Predictor: age, stress
    - Weight: weight | Cluster: psu

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    stress: HR = 1.367, p < .001 | weighted Cox + cluster-robust SE (weight, cluster=psu)
    time = followup, event = died ~ age + stress

>> COMMENTARY (narration):
    We ran the survival analysis under a complex sample: the trees came from a weighted and clustered design.
    Survey-PHREG builds a weighted Cox regression and gives cluster-robust standard errors; ignoring the design would
    lead to biased hazard ratios. The result is significant: the stress score raises the death hazard to 1.37 times
    (p < .001). What matters is that the estimate is population-generalizable -- because the design weights and
    clustering are handled correctly, the result represents the whole population, not the sample. In public-health
    cohorts and stratified RCT analyses, such design-aware survival models are essential; otherwise the standard
    errors mislead and results are biased.

====================================================================================

#92  Interval-Censored Survival
    file: 92_interval_censored_interval_takip.xlsx
  >> SCENARIO (narration):
    We followed trees with periodic inspections; the exact time of death is unknown -- we only know it occurred between two
    inspections, in [lower_bound_year, upper_bound_year]. This interval censoring misleads classical Kaplan-Meier. Interval-
    censored survival (Turnbull) handles this uncertainty correctly to estimate the survival curve: what is the median lifespan,
    and does it change by treatment (treatment)?
  >> VARIABLE SELECTION:
    - Lower bound (L): lower_bound_year
    - Upper bound (R): upper_bound_year

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Median survival = 10.0 | Turnbull (interval-censored Kaplan-Meier)
    event time within [lower_bound_year, upper_bound_year]

>> COMMENTARY (narration):
    We followed trees with periodic inspections; the exact time of death is unknown -- we only know it occurred
    within an interval between two inspections. This interval censoring misleads classical Kaplan-Meier, which assumes
    an exact event time. Interval-censored survival (the Turnbull algorithm) handled this uncertainty correctly to
    estimate the survival curve: a median lifespan of 10 time units. In the real world, data are rarely perfectly
    timed -- trees are checked at intervals, not every day. The value of this method is being able to model this
    realistic uncertainty without ignoring it, on a sound statistical basis. Forcing it into classical KM would have
    systematically mis-estimated the median lifespan.

====================================================================================

#93  Frailty Cox
    file: 93_frailty_cox_herd_based.xlsx
  >> SCENARIO (narration):
    Trees belong to stand/group clusters (herd_id); trees within the same group share conditions and therefore have similar death
    risk. Ignoring this within-group correlation would be biased. Frailty Cox adds a random "frailty" effect to each group: it
    correctly estimates the effect of stress on death risk while also modeling the hidden heterogeneity at the group level.
  >> VARIABLE SELECTION:
    - Time: followup
    - Event: died
    - Predictor: age, stress
    - Frailty (group): herd_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    stress: HR = 1.202, p < .001 | cluster-based frailty (herd_id), cluster-robust
    time = followup, event = died ~ age + stress

>> COMMENTARY (narration):
    Trees belonged to stand/group clusters; trees within the same group share conditions and thus have similar death
    risk. Ignoring this within-group correlation would be biased. Frailty Cox modeled the hidden group-level
    heterogeneity to correctly estimate the effect of stress on death risk: HR = 1.20, p < .001. So even after
    accounting for group differences, stress significantly raises death risk. This is the correct approach for
    hierarchical/clustered survival data -- trees in the same stand, animals in the same herd, patients in the same
    hospital; it includes the shared group risk in the model without breaking the independence assumption.

====================================================================================

#94  Time Series Analysis
    file: 94_time_series_monthly_insect.xlsx
  >> SCENARIO (narration):
    We examine the course of monthly insect counts (insect_count) over time (date). Time series analysis decomposes the series
    into trend, seasonal and irregular components and shows the autocorrelation structure: does the insect population follow a
    regular seasonal cycle within the year, is there a rising trend? Temperature/precipitation can be examined as extra
    explanatory variables.
  >> VARIABLE SELECTION:
    - Date: date
    - Series (value): insect_count

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    60 monthly observations | ADF p = 0.88 (non-stationary) | Trend: none | Seasonality: DETECTED
    insect_count time series

>> COMMENTARY (narration):
    We examined the course of monthly insect counts over time. Time series analysis decomposes the series into its
    components and tests the autocorrelation structure. Two important findings: first, the series is non-stationary
    (ADF p = 0.88) -- its statistical properties change over time, and differencing may be needed before modeling.
    Second and most interesting, seasonality was detected. So the insect population follows a regular cycle within the
    year -- rising in certain months, falling in others. There is no clear long-term trend. In practice this says the
    insect monitoring and control schedule should be planned around the seasonal peak: positioning intervention just
    before the peak months is the most effective strategy.

====================================================================================

#95  STL Decomposition (Seasonal Decomposition)
    file: 95_stl_pm25_mevsimsel_ayrisma.xlsx
  >> SCENARIO (narration):
    We want to decompose the monthly PM25 series (date, PM25) into its components: trend, seasonal component and residual. STL
    (LOESS-based) decomposition clearly shows how seasonality repeats within the year and the long-term trend: does air pollution
    peak in winter, and is there a rising tendency across the years?
  >> VARIABLE SELECTION:
    - Date: date
    - Series (value): PM25

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    STL (period = 7) | components: trend + seasonal + residual
    PM25 series decomposed

>> COMMENTARY (narration):
    We decomposed the PM25 (air pollution) series into its components: trend, seasonal component and residual. STL
    (LOESS-based) decomposition clearly shows how seasonality repeats within the cycle and the long-term trend. With
    period 7, a weekly/periodic seasonal pattern was captured. The practical value of this decomposition is large:
    looking at the raw series, saying "is pollution rising" can be misleading -- because seasonal fluctuation hides
    the trend. STL separates the two to reveal the genuine tendency. We can thus clearly answer "is air pollution
    seasonal, or is there a genuinely rising trend"; policy and warning systems are built on this distinction.

====================================================================================

#96  ARIMA / SARIMA (Forecasting)
    file: 96_arima_monthly_inflation.xlsx
  >> SCENARIO (narration):
    We want to forecast the monthly inflation series (inflation_monthly_pct) into the future. ARIMA exploits the series' own past
    values and error structure (auto-regressive + moving average); if seasonality is present, the SARIMA extension is used. The
    model produces point forecasts and confidence intervals for the coming months: where is the inflation trend heading?
  >> VARIABLE SELECTION:
    - Date: date
    - Series (value): inflation_monthly_pct

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    ARIMA   AIC = 197.3 | point forecasts + confidence intervals for coming months
    inflation_monthly_pct series

>> COMMENTARY (narration):
    We wanted to forecast the monthly inflation series into the future. ARIMA exploits the series' own past values
    and error structure (auto-regressive + moving average). The model fit the observed data well (AIC = 197.3) and
    produced point forecasts with confidence intervals for the coming months. ARIMA's strength is capturing the past
    pattern and projecting it forward; it is used directly in economic decision-making, budget planning and inventory
    management. Important note: the forecast's confidence interval widens as the horizon lengthens -- the near future
    is more certain, the far future more uncertain; a fundamental fact to keep in mind in any forecast.

====================================================================================

#97  Exponential Smoothing (Holt-Winters)
    file: 97_ets_monthly_visitor.xlsx
  >> SCENARIO (narration):
    We forecast the monthly number of visitors (visitor); the data have both trend and seasonality. Exponential smoothing
    (Holt-Winters) gives exponentially decaying weight to past observations and updates the level, trend and seasonal components;
    it is simpler than ARIMA but very effective in seasonal forecasting: what will the next season's visitor count be?
  >> VARIABLE SELECTION:
    - Date: date
    - Series (value): visitor

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Holt-Winters (seasonal additive, period = 12)   AIC = 1689 | next-season forecast ~ 77,099 visitors
    visitor series

>> COMMENTARY (narration):
    We forecast the monthly number of visitors; the data had both trend and seasonality. Exponential smoothing
    (Holt-Winters) gives exponentially decaying weight to past observations and updates the level, trend and seasonal
    components -- weighting the recent past more. It is simpler than ARIMA but very effective in seasonal forecasting.
    The model predicts about 77 thousand visitors for the next period. The practical value is large: staffing,
    parking and infrastructure planning of forest recreation areas rest on seasonal visitor forecasts. Knowing the
    peak season in advance lets resources be positioned at the right time and place -- neither short nor excess.

====================================================================================

#98  Mann-Kendall Trend Test
    file: 98_mann_kendall_temperature_trend.xlsx
  >> SCENARIO (narration):
    We test whether the annual mean temperature (annual_mean_temperature) shows a significant trend over the years (year).
    Mann-Kendall is a non-parametric trend test: without any distributional assumption it indicates whether the series has a
    monotonic (increasing/decreasing) tendency and, with Sen's slope, the rate of change: is the climate-warming signal
    statistically significant?
  >> VARIABLE SELECTION:
    - Time: year
    - Series (value): annual_mean_temperature

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    p < .001 ***   Trend = ↑ INCREASING   Sen's slope = 0.027 units/year
    annual_mean_temperature ~ year

>> COMMENTARY (narration):
    We tested whether the annual mean temperature shows a significant trend over the years. Mann-Kendall is a
    non-parametric trend test: without any distributional assumption it checks for a monotonic tendency. The result is
    clear and important: p below one in a thousand, with an upward trend -- a statistically significant INCREASE in
    temperature. Sen's slope gives the rate of this rise robustly: about 0.027 units per year. This is a classic
    climate-warming signal. Mann-Kendall is preferred because climate/environmental series are usually non-normal and
    contain outliers; this test is resistant to them. In practice such trend analyses form the basis of climate-change
    monitoring and long-term forest-management strategies.

====================================================================================

#99  Anomaly Detection (Outliers)
    file: 99_anomali_detection_outlier_inventory.xlsx
  >> SCENARIO (narration):
    We want to automatically detect abnormal/erroneous records in the inventory measurements (dbh, height, biomass). Anomaly
    detection flags individuals that deviate markedly from the others in the multivariate space (e.g. a measurement error or an
    extraordinary tree): which records should be reviewed for data cleaning -- a genuine outlier or an error?
  >> VARIABLE SELECTION:
    - Variables: dbh, height, biomass

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Anomaly rate = 8.5% (IQR method) | lower/upper bound + anomaly count per variable
    dbh, height, biomass

>> COMMENTARY (narration):
    We wanted to automatically detect abnormal/erroneous records in the inventory measurements. Anomaly detection
    flags individuals that deviate markedly from the others in the multivariate space. With the IQR method, 8.5
    percent of observations were flagged as anomalies. What matters is distinguishing "error" from "genuine outlier":
    a measurement error (e.g. 80 m height) should be corrected, but an extraordinarily large monument tree is real and
    valuable. In practice this is the first step of data cleaning: the flagged records are reviewed manually. Automatic
    anomaly detection is the practical way to systematically catch erroneous records that the eye would miss in
    thousands of rows.

====================================================================================

#100  Variance Components (h2)
    file: 100_varcomp_3seviye_h2.xlsx
  >> SCENARIO (narration):
    In a three-level nested trial (upper_unit / lower_unit / measurement) we want to decompose from which level the total variance
    of a trait (value) originates. Variance-components analysis partitions the variance into levels; in forestry/breeding this is
    the basis for heritability (h2) estimation: how much of the variability in the trait is genetic/upper-group, and how much is
    environmental/measurement?
  >> VARIABLE SELECTION:
    - Dependent (continuous): value
    - Upper level: upper_unit
    - Lower level (nested): lower_unit

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Variance components (REML):  upper_unit 12.1%  |  lower_unit (nested) 40.9%  |  residual 47.0%
    dv = value

>> COMMENTARY (narration):
    In a nested trial we decomposed from which level the total variance of a trait originates. Variance-components
    analysis, unlike LMM, uses no fixed effects -- all factors are random, the goal being to partition variance into
    levels. The result is instructive: 12 percent of the variability comes from the upper level (e.g. population), 41
    percent from the lower level (e.g. family, nested within population), and 47 percent from the residual
    (environment/measurement). In forestry and breeding this is the basis of heritability (h2) estimation: how much of
    the variability in a trait is genetic/structural and how much environmental? The lower level (family) contributing
    more than the upper (population) directly shows which level to focus on in selection/breeding work.

====================================================================================

#101  Bayesian t-Test
    file: 101_bayesian_t_test_new_old.xlsx
  >> SCENARIO (narration):
    We test the difference between two applications, new and old (group), on a measurement (value) from a Bayesian framework. The
    classical t-test gives only a p-value; the Bayesian t-test, with the Bayes factor (BF10), answers "by how many times do the
    data support one hypothesis over the other" and shows the posterior distribution of the effect size: how strong is the evidence
    for a difference?
  >> VARIABLE SELECTION:
    - Dependent (continuous): value
    - Grouping (2 categories): group

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    BF10 = 3.4 x 10^46 (decisive evidence, H1)   Cohen's d = 3.94   (value ~ group: new/old)

>> COMMENTARY (narration):
    We tested the difference between a new and an old application on a measurement from a Bayesian framework. The
    classical t-test gives only a p-value; the Bayesian t-test, with the Bayes factor, answers "by how many times do
    the data support one hypothesis over the other". The result is staggering: BF10 = 3.4 x 10 to the 46. This means
    the data support "a difference exists" over "no difference" by an astronomical ratio -- far beyond "decisive
    evidence". Cohen's d = 3.94 means the effect is enormous too. The beauty of the Bayes factor is doing what the
    p-value cannot: it can also measure evidence FOR H0 and gives the STRENGTH of evidence on a continuous scale. Here
    the difference is so clear that we can say the new application is decisively superior to the old.

====================================================================================

#102  Bayesian Correlation
    file: 102_bayesian_correlation_BF10.xlsx
  >> SCENARIO (narration):
    We examine the relationship between two continuous variables (X_variable, Y_variable) in a Bayesian way. Classical correlation
    gives an r and a p; Bayesian correlation provides the posterior distribution of r and the Bayes factor (BF10): how strong is the
    evidence for "a relationship exists" against "no relationship", and what is the credible interval for r?
  >> VARIABLE SELECTION:
    - Variable 1: X_variable
    - Variable 2: Y_variable

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    r = 0.780 (very strong)   BF10 = 4.5 x 10^14 (decisive: relationship exists)   (X_variable, Y_variable)

>> COMMENTARY (narration):
    We examined the relationship between two continuous variables in a Bayesian way. Classical correlation gives an r
    and a p; Bayesian correlation provides the posterior distribution of r and the Bayes factor. The result is very
    strong: r = 0.78 and BF10 = 4.5 x 10 to the 14. So "a relationship exists" is supported over "no relationship" by
    a fourteen-zero ratio -- decisive evidence. The advantage of the Bayesian approach is answering not just "is it
    significant" but "how strong is the evidence" and "what is a reasonable range for r" directly. In a very strong
    relationship like this, the Bayes factor makes the result indisputable: there is a real, strong connection between
    the two variables.

====================================================================================

#103  Bayesian ANOVA
    file: 103_bayesian_anova_2yonlu.xlsx
  >> SCENARIO (narration):
    We examine the effect of two factors (factor1, factor2) on an outcome (value) with a Bayesian two-way ANOVA. Classical ANOVA
    gives F and p; Bayesian ANOVA compares which model (main effects, interaction) best explains the data using Bayes factors: is
    there evidence for an interaction, and which effect is most strongly supported?
  >> VARIABLE SELECTION:
    - Dependent (continuous): value
    - Factor 1: factor1
    - Factor 2: factor2

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Bayesian ANOVA | factor(s): factor1 | model comparison (Bayes factors)
    value ~ factor1 (+ factor2)

>> COMMENTARY (narration):
    We examined the effect of factor(s) on an outcome with a Bayesian ANOVA. Classical ANOVA gives F and p; Bayesian
    ANOVA compares which model -- main effects only, or also the interaction -- best explains the data using Bayes
    factors. This gives a direct answer to "is there evidence for an interaction"; in the classical approach a
    "non-significant" interaction does NOT prove its absence, whereas the Bayes factor can measure evidence for H0 too.
    The practical value of Bayesian ANOVA is expressing model selection in the language of probability: being able to
    say "the data most support this model" is more informative for decision-making than classical p-value thresholds.

====================================================================================

#104  Hierarchical Bayesian (LMM)
    file: 104_hierarchical_bayesian_LMM.xlsx
  >> SCENARIO (narration):
    In data nested within groups (group_id) we estimate the effect of a covariate (X_covariate) on an outcome (Y_response) with a
    hierarchical Bayesian model. The hierarchical Bayesian approach estimates group-level effects with partial pooling: groups with
    little data are pulled toward the overall mean, yielding more stable estimates that express uncertainty explicitly.
  >> VARIABLE SELECTION:
    - Dependent (continuous): Y_response
    - Covariate (fixed effect): X_covariate
    - Grouping (random): group_id

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    σ²_u (group random) = 2.90 (8.6%) | X_covariate fixed effect with 95% credible interval
    Y_response ~ X_covariate + (1 | group_id)

>> COMMENTARY (narration):
    In data nested within groups we estimated the effect of a covariate on an outcome with a hierarchical Bayesian
    model. The hierarchical Bayesian approach estimates group-level effects with "partial pooling": groups with little
    data are pulled toward the overall mean, yielding more stable estimates that express uncertainty explicitly. The
    group random variance is 8.6 percent of the total variability -- so there is a notable but not dominant difference
    among groups. The covariate's effect was reported with a posterior credible interval. The power of this approach
    shows especially in low-data/many-group situations: it "intelligently" corrects individual group estimates,
    prevents overfitting, and conveys uncertainty honestly.

====================================================================================

#105  Spatial Lag (SAR)
    file: 105_spatial_sar_spatial.xlsx
  >> SCENARIO (narration):
    Across spatial units (lat, lon) we model an outcome (Y_value) with X1, X2; but neighboring units influence one another (spatial
    spillover). Classical regression ignores this dependence. The spatial-lag model (SAR) brings the neighbors' Y into the model to
    capture the spillover: does a high value in one unit raise its neighbors too (spatial spillover effect)?
  >> VARIABLE SELECTION:
    - Dependent (outcome): Y_value
    - Predictors: X1, X2
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    ρ (rho) = 0.336   AIC = 1051.6 | k-NN W matrix (Y_value ~ X1 + X2, lat/lon)

>> COMMENTARY (narration):
    We modeled an outcome with X1, X2 across spatial units; but neighboring units influence one another (spatial
    spillover). Classical regression ignores this dependence. The spatial-lag model (SAR) captured the spillover by
    bringing the neighbors' Y into the model: the spatial autoregressive parameter rho = 0.34, positive and notable.
    This means a high value in one unit raises its neighbors too -- a spatial "spillover" effect. In practice this is
    very important: e.g. if a high-yield parcel also affects neighboring parcels' yield, interventions should be planned
    not one by one but in spatial clusters. SAR prevents the biased estimates that ignoring spatial dependence would cause.

====================================================================================

#106  Spatial Error (SEM)
    file: 106_spatial_error_residual.xlsx
  >> SCENARIO (narration):
    Again across spatial units we model Y_value with X1, X2; but this time the spatial dependence is in the residuals (from
    unobserved shared environmental factors). The spatial-error model (SEM) models the spatial structure of the error term to keep
    the coefficients from becoming biased/inflated: what are the X effects once they are purified from the hidden spatial residual
    arising from the neighborhood?
  >> VARIABLE SELECTION:
    - Dependent (outcome): Y_value
    - Predictors: X1, X2
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    λ (lambda) = 0.590   z = 6.38, p < .001   Pseudo R² = 0.84   AIC = 1037.4
    Y_value ~ X1 + X2 (lat/lon)

>> COMMENTARY (narration):
    Again across spatial units we modeled Y with X1, X2; but this time the spatial dependence is in the residuals (from
    unobserved shared environmental factors). The spatial-error model (SEM) models the spatial structure of the error
    term to keep the coefficients from becoming biased/inflated. The result is significant: lambda = 0.59, z = 6.38, p <
    .001 -- a strong spatial autocorrelation in the residuals. This supports the "omitted spatial variables" hypothesis:
    something spatially structured that we did not include is affecting the residuals. With pseudo R² = 0.84 the model is
    still strong. Choosing between SAR and SEM means asking whether the spillover is in the outcome (SAR) or in the error
    (SEM); here AIC finds SEM slightly better.

====================================================================================

#107  GWR (Geographically Weighted Regression)
    file: 107_gwr_local.xlsx
  >> SCENARIO (narration):
    We think the effect of X1, X2 on Y_value is not the same everywhere but varies by location. Classical regression gives a single
    global coefficient; GWR estimates local coefficients for each location: a variable's effect may be strong in one region and
    weak/reversed in another. It makes spatial heterogeneity visible on the map.
  >> VARIABLE SELECTION:
    - Dependent (outcome): Y_value
    - Predictors: X1, X2
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    R² = 0.887   Adj. R² = 0.868   AICc = 1129.9 | local coefficients per location
    Y_value ~ X1 + X2 (lat/lon)

>> COMMENTARY (narration):
    We thought the effect of X1, X2 on Y is not the same everywhere but varies by location. Classical regression gives a
    single global coefficient; GWR estimates local coefficients for each location. The model fit very well (R² = 0.89).
    GWR's key output is the VARIATION of the coefficients across the map: a variable's effect may be strong in one
    region, weak or reversed in another. This makes "spatial heterogeneity" visible. In practice this is very valuable:
    e.g. if precipitation's effect on yield is strong in an arid region and weak in a humid one, a single global
    coefficient hides it -- GWR reveals it. Far more informative than global models for region-specific management decisions.

====================================================================================

#108  Nested LMM
    file: 108_nested_lmm_R_P_F.xlsx
  >> SCENARIO (narration):
    In a replicated hierarchical breeding trial -- block (block), upper group (upper_group), lower group (lower_group_no) -- we
    estimate the variance components and effects of a trait (value) in a single model. The nested LMM, in this structure where the
    lower group is nested within the upper group, gives each level's contribution with correct F tests: is the difference among
    populations/families significant, is there a genotype-environment interaction?
  >> VARIABLE SELECTION:
    - Dependent (continuous): value
    - Replication (block): block
    - Upper group: upper_group
    - Lower group (nested): lower_group_no

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    upper_group (P): F(3) = 4.92, p = 0.010 * | model: value ~ R + P + R*P + F(P) + R*F(P)
    block (R), upper_group (P), lower_group_no (F nested in P)

>> COMMENTARY (narration):
    In a replicated hierarchical breeding trial -- block, upper group (population), lower group (family, nested within
    population) -- we estimated the variance components and effects of a trait in a single model. The nested LMM, in
    this nested structure, gives each level's contribution with correct F tests (each effect uses the appropriate error
    term -- something ordinary ANOVA cannot do). The result: the upper group (population) is significant, F(3) = 4.92, p
    = 0.010. So there is a genetic/structural difference among populations. This is the fundamental question of forest
    breeding: at which level should selection be done? If the population difference is significant, choosing the right
    population is the priority; if the family difference dominates, choosing families within a population is. The nested
    LMM provides the statistical basis for this decision.

====================================================================================

#109  Crossed LMM
    file: 109_crossed_lmm_A_B.xlsx
  >> SCENARIO (narration):
    In a design where two random factors (factor_a, factor_b) are crossed -- each A level pairs with each B level -- we model the
    outcome (value) (e.g. every genotype tested at every location). The crossed LMM takes A and B as crossed random effects; it thus
    models the variability from both genotype and location at once and gives purified estimates.
  >> VARIABLE SELECTION:
    - Dependent (continuous): value
    - Crossed random effects: factor_a, factor_b

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    factor_a (random): F(3) = 0.38, p = 0.767 ns | value ~ factor_a + factor_b (crossed random)

>> COMMENTARY (narration):
    In a design where two random factors are crossed -- each A level pairs with each B level -- we modeled the outcome
    (e.g. every genotype tested at every location). The crossed LMM takes A and B as crossed random effects, modeling
    the variability from both genotype and location at once. The result: factor A is not significant (F(3) = 0.38, p =
    0.77) -- so there is no notable difference among its levels. Non-significance is also a finding: factor A's
    contribution to the outcome is negligible. Crossed designs differ from nested ones: here A and B are independently
    crossed (every combination present), whereas in nested the sub-factor is embedded in the upper. The correct model
    depends on the true structure of the design.

====================================================================================

#110  KDE (Kernel Density Map)
    file: 110_KDE_field_distribution_density.xlsx
  >> SCENARIO (narration):
    We want to see the spatial density of field points (lat, lon) and yield (yield_ton_ha) as a continuous surface. Kernel density
    estimation (KDE) turns point data into a smooth density map: where does production concentrate, which areas are sparse? Hot and
    cold regions are read at a glance.
  >> VARIABLE SELECTION:
    - Coordinates: lat, lon
    - Weight (value): yield_ton_ha

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Map: kernel density (KDE) surface | points: lat/lon, weight: yield_ton_ha (visual output)

>> COMMENTARY (narration):
    We saw the spatial density of field points and yield as a continuous surface. Kernel density estimation (KDE) turns
    scattered point data into a smooth heat map: where production concentrates and which areas are sparse is read at a
    glance. This is a visual exploration tool -- instead of a single number, it shows the WHOLE spatial pattern. In
    precision agriculture this is very valuable: high-yield "hot" zones and low-yield "cold" zones can be distinguished
    and interventions (fertilization, irrigation) directed in a targeted way. KDE is the visual answer to "where does
    yield concentrate".

====================================================================================

#111  Hexbin (Density Map)
    file: 111_Hexbin_plot_distribution.xlsx
  >> SCENARIO (narration):
    We aggregate many field points (lat, lon) into hexagonal bins to show density. The hexbin map prevents point piles from
    overlapping in the display and shows, by color, how many observations / what density fall in each hexagon: the dense and sparse
    regions of the distribution become clear.
  >> VARIABLE SELECTION:
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Map: hexagonal-bin (hexbin) density map | points: lat/lon (visual output)

>> COMMENTARY (narration):
    We aggregated many field points into hexagonal bins to show density. The hexbin map solves the overlap problem: when
    thousands of points hide one another, hexbin shows by color how many observations fall in each hexagon. Unlike KDE,
    it uses discrete bins instead of a continuous surface -- which lets density be read numerically and clearly. In
    practice it is the cleanest way to visualize large point datasets (e.g. all parcel locations): the dense and sparse
    regions, the structure of the distribution, become clear at a glance. Whenever a scatter plot becomes too crowded,
    hexbin is preferred.

====================================================================================

#112  Moran's I (Spatial Autocorrelation)
    file: 112_Morans_I_yield_autocorrelation.xlsx
  >> SCENARIO (narration):
    Is field yield (yield_ton_ha) distributed spatially at random, or are similar values clustered near one another? Moran's I
    measures spatial autocorrelation in a single index: a positive I indicates that high-yield fields are adjacent to high-yield
    neighbors (clustering). Is there a spatial pattern in yield, and is it significant?
  >> VARIABLE SELECTION:
    - Value: yield_ton_ha
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Moran's I = 0.405   p = 0.001 ** (positive spatial autocorrelation)   yield_ton_ha (lat/lon, k-NN W)

>> COMMENTARY (narration):
    We tested whether field yield is distributed spatially at random, or whether similar values cluster near one
    another. Moran's I measures spatial autocorrelation in a single index (-1 to +1; 0 = random). The result: I = 0.405,
    p = 0.001 -- positive and significant. The meaning is clear: high-yield fields are adjacent to high-yield neighbors,
    and low-yield ones cluster too; yield is spatially CLUSTERED, not randomly distributed. This is a very important
    finding: if a non-random spatial pattern exists, classical statistics (assuming independence) become misleading and
    spatial models (SAR/SEM) are needed. In practice this clustering in yield points to an underlying spatial gradient
    (soil, water, microclimate) and opens the door to targeted land management.

====================================================================================

#113  Getis-Ord G* (Hotspot Analysis)
    file: 113_Getis_Ord_high_yield_hotspot.xlsx
  >> SCENARIO (narration):
    We want to find where yield (yield_ton_ha) forms statistically significant "hotspots" (high-high clusters) and "coldspots". For
    each location, Getis-Ord G* compares the local sum of surrounding values: is this field at the center of a high-yield cluster?
    Which parcels are priorities for targeted management?
  >> VARIABLE SELECTION:
    - Value: yield_ton_ha
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Map: local G* z-scores (hot/cold spots) | yield_ton_ha (lat/lon)
    (global clustering confirmed by Moran's I = 0.41, p = .001)

>> COMMENTARY (narration):
    We examined where yield forms statistically significant "hotspots" (high-high clusters) and "coldspots". Moran's I
    tells us there is GLOBAL clustering; Getis-Ord G* LOCALIZES it: for each location it compares the local sum of
    surrounding values to tell whether that point is part of a hot or cold cluster, via a z-score. The map shows
    high-yield hotspots in red and low-yield coldspots in blue. The practical value is direct: hotspots can be studied
    for "why so productive" and good practices propagated; coldspots are areas needing priority intervention (soil
    improvement, drainage). The Moran's + Getis-Ord pair answers "is there clustering" and "where" together.

====================================================================================

#114  DBSCAN (Spatial Field Clusters)
    file: 114_DBSCAN_field_kumeleri.xlsx
  >> SCENARIO (narration):
    We want to split field/parcel locations (lat, lon) into natural clusters by geographic proximity. DBSCAN performs density-based
    clustering: nearby parcels form a cluster, isolated parcels are separated as "noise". There is no need to specify the number of
    clusters in advance; it directly reveals irregularly shaped agricultural-area clusters and outlier locations.
  >> VARIABLE SELECTION:
    - Coordinates: lat, lon

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    DBSCAN: 4 clusters + 11 noise points (n = 115) | lat/lon, density-based

>> COMMENTARY (narration):
    We split field/parcel locations into natural clusters by geographic proximity. DBSCAN performs density-based
    clustering: nearby parcels form a cluster, isolated parcels are separated as "noise". In this data, 4 natural field
    clusters and 11 isolated (noise) parcels were found. It differs from k-means in two key ways: no need to specify the
    number of clusters in advance, and it can capture irregularly shaped clusters -- plus it explicitly flags outlier
    locations (noise). In practice this is ideal for dividing scattered agricultural land into management units: each
    cluster can be treated as an operational block, and isolated parcels planned separately. Directly usable for spatial
    organization and logistics planning.

====================================================================================

#115  Mixed-Design (Split-Plot) ANOVA
    file: 115_mixed_anova_crop_yield_t_ha.xlsx
  >> SCENARIO (narration):
    We follow 40 plots measured at three season levels (season1, season2, season3); each belongs to one of two fertilizer groups (organic / mineral). A mixed (split-plot) design tests the between-subjects
    main effect, the within-subjects main effect and their interaction on crop yield.
  >> VARIABLE SELECTION:
    - Dependent variable: crop_yield_t_ha
    - Subject ID: plot_id
    - Between-subjects factor: fertilizer
    - Within-subjects factor: season

  ----- RESULT & COMMENTARY -----
>> RESULT (screen):
    Between (fertilizer): F(1,38) = 5.76  p = 0.021  np2 = 0.132
    Within (season):  F(2,76) = 148.64  p < .001  np2 = 0.796
    Interaction:         F(2,76) = 11.09  p < .001  np2 = 0.226
    Mauchly W = 0.781  p = 0.009   n_subjects = 40   n_obs = 120

>> COMMENTARY (narration):
    In a fertilizer x season mixed design we analyzed crop yield for 40 plots (120 observations). The interaction is significant (F(2,76) = 11.09, p < .001, np2 = 0.226) *** -- the two groups' change across season differs in magnitude. The between-subjects main effect (organic vs mineral) is F = 5.76, p = 0.021; the within-subjects main effect (season1/season2/season3) is F = 148.64, p < .001. Mauchly's test p = 0.009, so sphericity is violated, so the Greenhouse-Geisser corrected within p is read. Read the interaction first: when it is significant the group effect must be interpreted separately at each season level. In agriculture, the mixed design is the standard analysis for comparing fertilizer regimes on yield across seasons.

====================================================================================
