Reconciling missing stems with apparent recruits in multi-census tree inventory data

2026-09-18

When a tree inventory census is completed in a long-term plot, looking back at the data later the situation may arise where some stems from the previous census have seemingly disappeared. This could happen naturally through stem mortality and subsequent decomposition. At the same time, there may be many new stems that have apparently recruited during the most recent census. Some of these apparent recruits may actually be the missing stems from the previous census, for example if the tags have dropped off. Matching the missing stems with the apparent recruits can be a laborious and complex puzzle.

I wrote a function that compares missing trees with apparent recruit records to evaluate the likelihood of an apparent recruit actually being a missing stem that was recorded in the previous census. The function calculates a probabilistic “trust score” based on spatial proximity, diameter growth trajectories, and taxonomic IDs. The function identifies complex matching topologies: 1-to-1, 1-to-many, many-to-1, and many-to-many.

Here is the function:

#' Match missing trees with apparent recruits in re-census data
#' 
#' Compares trees from a previous census that are missing in the current
#' census, with trees that apparently recruited in the current census to
#' determine whether any of the apparent recruits may actually be the missing
#' trees from the previous census. Evaluates spatial proximity, diameter growth
#' trajectories, and taxonomy to calculate a "trust score" for potential
#' matches. Handles complex topologies (1-to-1, 1-to-many, many-to-1) with
#' optional common sense filters to exclude small diameter true recruits from
#' the matching process.
#'
#' @param miss data frame containing trees that were identified as
#'     "missing" in the current census. Measurement values should be from the
#'     previous census. 
#' @param miss_id column name in `miss` containing tree IDs 
#' @param miss_x column name in `miss` containing X spatial coordinate values
#' @param miss_y column name in `miss` containing Y spatial coordinate values
#' @param miss_diam column name in `miss` containing diameter measurement values
#' @param miss_genus column name in `miss` containing genus names
#' @param miss_species column name in `miss` containing binomial species names
#' 
#' @param rec data frame containing trees that were identified as "recruits" in
#'     the current census. Measurement values should be from the current
#'     census.
#' @param rec_id column name in `rec` containing tree IDs 
#' @param rec_x column name in `rec` containing X spatial coordinate values
#' @param rec_y column name in `rec` containing Y spatial coordinate values
#' @param rec_diam column name in `rec` containing diameter measurement values
#' @param rec_genus column name in `rec` containing genus names
#' @param rec_species column name in `rec` containing binomial species names
#' 
#' @param max_dist maximum spatial distance threshold to consider a match.
#' @param sd_dist spatial decay parameter representing trust in coordinate
#'     precision (expected spatial coordinate uncertainty). Lower values assume
#'     coordinates are very accurate, such that only missing stems very close
#'     to recruits are considered as potential matches.
#' 
#' @param inc_diam expected average diameter growth increment between the two
#'     censuses.
#' @param sd_inc_diam tolerance parameter for diameter growth variation.
#'     Controls how forgiving the model is of individual trees growing faster
#'     or slower than expected
#' @param max_neg_diam maximum allowable negative diameter difference (i.e.,
#'     shrinkage) before applying a heavy penalty for impossible shrinkage.
#' @param max_rec_diam Maximum expected diameter of true recruits. Only trees
#'     larger than this are considered as potential missing individuals. Set to
#'     zero to consider all apparent recruits as potential missing trees.
#' 
#' @param min_trust minimum trust score required to consider a potential match.
#' 
#' @return A list containing `pairs` (all valid scored edges), `clusters` (group summaries 
#'         with match topologies), and `true_recruits` (unmatched true recruits).
#' 
recruitMissingMatch <- function(
  miss, miss_id, miss_x, miss_y, miss_diam, miss_genus, miss_species,
  rec, rec_id, rec_x, rec_y, rec_diam, rec_genus, rec_species, 
  max_dist = 5.0, sd_dist = 1.5, 
  inc_diam = 1.0, sd_inc_diam = 2.0, max_neg_diam = 2.0, max_rec_diam = 5.0, 
  min_trust = 0.02) {
  
  # Edge case handling 
  if (nrow(miss) == 0) { 
    stop("Data frame 'miss' contains no candidate missing trees")
  } 

  if (nrow(rec) == 0) { 
    stop("Data frame 'rec' contains no candidate missing trees")
  } 
  
  # Filter to potentially bad recruits based on large diameter
  is_potential_bad_rec <- rec[[rec_diam]] > max_rec_diam
  potential_bad_rec <- rec[is_potential_bad_rec, , drop = FALSE]
  true_recruits <- rec[!is_potential_bad_rec, , drop = FALSE] 
  
  if (nrow(potential_bad_rec) == 0) {
    message("After filtering out true recruits by 'max_rec_diam', no potential bad recruits remain")
    return(list(
      pairs = data.frame(), 
      clusters = data.frame(), 
      true_recruits = true_recruits
    ))
  }

  # Extract vectors using column name arguments 
  m_id   <- miss[[miss_id]]
  m_x    <- miss[[miss_x]]
  m_y    <- miss[[miss_y]]
  m_diam <- miss[[miss_diam]]
  m_gen  <- miss[[miss_genus]]
  m_sp   <- miss[[miss_species]]
  
  r_id   <- potential_bad_rec[[rec_id]]
  r_x    <- potential_bad_rec[[rec_x]]
  r_y    <- potential_bad_rec[[rec_y]]
  r_diam <- potential_bad_rec[[rec_diam]]
  r_gen  <- potential_bad_rec[[rec_genus]]
  r_sp   <- potential_bad_rec[[rec_species]]
  
  # Pairwise feature extraction and scoring
  pairs <- expand.grid(
    missing_idx = seq_len(nrow(miss)),
    recruit_idx = seq_len(nrow(potential_bad_rec))
  )
  
  mi <- pairs$missing_idx
  ri <- pairs$recruit_idx
  
  # Spatial component (Gaussian decay based on coordinate uncertainty)
  dist_val <- sqrt((m_x[mi] - r_x[ri])^2 + (m_y[mi] - r_y[ri])^2)
  score_dist <- exp(-(dist_val^2) / (2 * sd_dist^2))
  score_dist[dist_val > max_dist] <- 0
  
  # Diameter component (segregate positive growth vs negative shrinkage)
  diam_diff <- r_diam[ri] - m_diam[mi]
  score_diam <- numeric(length(diam_diff))
  
  pos_idx <- diam_diff >= 0
  score_diam[pos_idx] <- exp(-((diam_diff[pos_idx] - inc_diam)^2) / (2 * sd_inc_diam^2))
  
  neg_idx <- diam_diff < 0
  score_diam[neg_idx & diam_diff >= -max_neg_diam] <- 0.85
  score_diam[neg_idx & diam_diff < -max_neg_diam]  <- 0.05
  
  # Taxonomic component (genus/species matching)
  score_taxa <- ifelse(m_gen[mi] == r_gen[ri] & m_sp[mi] == r_sp[ri], 1.0, 
    ifelse(m_gen[mi] == r_gen[ri], 0.7, 0.1))
  
  # Composite trust score
  trust_score <- score_dist * score_diam * score_taxa
  
  # Compile pairwise edge data frame
  edge_df <- data.frame(
    missing_id    = m_id[mi],
    recruit_id    = r_id[ri],
    distance      = round(dist_val, 2),
    diam_t1       = m_diam[mi],
    diam_t2       = r_diam[ri],
    diam_diff     = round(diam_diff, 2),
    genus_match   = (m_gen[mi] == r_gen[ri]),
    species_match = (m_sp[mi] == r_sp[ri]),
    score_dist    = round(score_dist, 4),
    score_diam    = round(score_diam, 4),
    score_taxa    = round(score_taxa, 4),
    trust_score   = round(trust_score, 4)
  )
  
  # Filter by minimum trust threshold
  edge_df <- edge_df[edge_df$trust_score >= min_trust, ]
  
  if (nrow(edge_df) == 0) {
    message("No potential matches found")
    return(list(
      pairs = data.frame(), 
      clusters = data.frame(), 
      true_recruits = rec
    ))
  }
  
  # Network component clustering 
  edge_df$m_node <- paste0("M_", edge_df$missing_id)
  edge_df$r_node <- paste0("R_", edge_df$recruit_id)
  
  nodes <- unique(c(edge_df$m_node, edge_df$r_node))
  comp  <- setNames(seq_along(nodes), nodes)
  
  changed <- TRUE
  while (changed) {
    changed <- FALSE
    for (i in seq_len(nrow(edge_df))) {
      c1 <- comp[edge_df$m_node[i]]
      c2 <- comp[edge_df$r_node[i]]
      if (c1 != c2) {
        min_c <- min(c1, c2)
        comp[comp == c1 | comp == c2] <- min_c
        changed <- TRUE
      }
    }
  }
  
  edge_df$cluster_id <- comp[edge_df$m_node]
  edge_df$m_node     <- NULL
  edge_df$r_node     <- NULL
  
  # Summarize cluster topologies
  cluster_list <- split(edge_df, edge_df$cluster_id)
  summary_rows <- vector("list", length(cluster_list))
  
  for (i in seq_along(cluster_list)) {
    sub_df   <- cluster_list[[i]]
    m_unique <- unique(sub_df$missing_id)
    r_unique <- unique(sub_df$recruit_id)
    n_m      <- length(m_unique)
    n_r      <- length(r_unique)
    
    if (n_m == 1 && n_r == 1) {
      match_type <- "1-to-1" 
    } else if (n_m == 1 && n_r > 1) {
      match_type <- "1-to-Many"
    } else if (n_m > 1 && n_r == 1) {
      match_type <- "Many-to-1"
    } else {
      match_type <- "Many-to-Many"
    }    

    summary_rows[[i]] <- data.frame(
      cluster_id = as.integer(names(cluster_list)[i]),
      match_type = match_type,
      missing_trees = paste(m_unique, collapse = ", "),
      recruit_trees = paste(r_unique, collapse = ", "),
      n_missing = n_m,
      n_recruits = n_r,
      mean_trust = round(mean(sub_df$trust_score), 4),
      min_trust = min(sub_df$trust_score),
      max_trust = max(sub_df$trust_score)
    )
  }
  
  clusters_summary <- do.call(rbind, summary_rows)
  rownames(clusters_summary) <- NULL
  clusters_summary <- clusters_summary[order(-clusters_summary$mean_trust), ]
  
  # Aggregate final true recruits
  matched_bad_rec_ids <- unique(edge_df$recruit_id)
  unmatched_bad_rec  <- potential_bad_rec[!potential_bad_rec[[rec_id]] %in% matched_bad_rec_ids, , drop = FALSE]
  final_true_recruits <- rbind(true_recruits, unmatched_bad_rec)
  
  return(list(
    pairs = edge_df[order(-edge_df$trust_score), ],
    clusters = clusters_summary,
    true_recruits = final_true_recruits
  ))
}

recruitMissingMatch() uses a set of configurable parameters to tweak the matching algorithm to the specific characteristics and measurement error rates of a given forest plot. The arguments are grouped into data inputs, column mappings, spatial bounds, diameter growth mechanics, and filtering thresholds:

Data inputs and column mappings:

The function requires two datasets and explicit column name mappings:

Spatial parameters:

These parameters govern the spatial decay function, which calculates the likelihood of a match based on the physical distance between coordinates.

Diameter growth parameters:

These parameters evaluate plausibility and field measurement error using a scoring system for positive growth versus negative shrinkage.

Filtering and network thresholds:

These parameters control which stems enter the network clustering phase and what constitutes a valid match.

Examples

The following examples demonstrate the function’s behavior across different data scenarios. The tibble::tribble() function is used to construct the test datasets.

missing_census <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,    ~species,
  "T_101", 10.2, 12.4,  18.0, "Quercus", "alba"
)

recruit_census <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,    ~species,
  "R_901", 10.5, 12.8,  19.2, "Quercus", "alba",
  "R_903", 10.8, 13.0,   3.2, "Quercus", "alba"
)

result <- recruitMissingMatch(
  miss = missing_census, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5, sd_dist = 1.5, inc_diam = 1, sd_inc_diam = 2,
  max_neg_diam = 2, max_rec_diam = 5, min_trust = 0.02)

result

Outcome: The function identifies one valid match (T_101 to R_901). The stem R_903 is excluded from matching and returned in the true_recruits data frame due to its small diameter.

Example 2: Handling 1-to-Many Ambiguity

When a missing stem is spatially proximate to multiple recruits, the function groups them for review.

missing_census_ambig <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,    ~species,
  "T_101", 50.0, 50.0,  20.0, "Acer",    "saccharum",
  "T_102", 20.5, 30.2,  12.5, "Quercus", "alba"
)

recruit_census_ambig <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,    ~species,
  "R_901", 50.4, 50.2,  21.2, "Acer",    "saccharum",
  "R_902", 51.2, 49.5,  23.0, "Acer",    "rubrum",
  "R_903", 20.6, 30.4,   3.2, "Quercus", "alba"
)

results_ambig <- recruitMissingMatch(
  miss = missing_census_ambig, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_ambig, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5, sd_dist = 1.5, inc_diam = 1.2, sd_inc_diam = 2,
  max_neg_diam = 2, max_rec_diam = 5, min_trust = 0.02)

results_ambig

Outcome: The function returns two potential matches for T_101. R_901 receives a higher trust score than R_902. R_903 is classified as a true recruit.

Example 3: Many-to-1 Consolidation

This scenario models two previously distinct stems that were consolidated into a single record in the current census. The inc_diam tolerance is increased to account for the artificially inflated expected growth.

missing_census_many <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 30.0, 30.0,  8.5,  "Pinus", "taeda",
  "T_102", 30.8, 30.2,  9.0,  "Pinus", "taeda"
)

recruit_census_many <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "R_901", 30.4, 30.1, 19.5,  "Pinus", "taeda"
)

results_many <- recruitMissingMatch(
  miss = missing_census_many, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_many, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 3, sd_dist = 1.5, inc_diam = 10, sd_inc_diam = 4,
  max_neg_diam = 2, max_rec_diam = 5, min_trust = 0.0001)

results_many

Outcome: Apparent recruit R_901 is matched with both T_101 and T_102, forming a “Many-to-1” cluster.

Example 4: Many-to-Many Complex Relationship

In dense clusters with identical species, multiple missing stems and recruits may overlap within the defined thresholds. The network clustering algorithm groups these into a single connected component.

missing_census_many_many <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 10.0, 10.0, 15.0,  "Fagus", "grandifolia",
  "T_102", 11.5, 10.5, 18.0,  "Fagus", "grandifolia"
)

recruit_census_many_many <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "R_901", 10.2, 10.2, 16.2,  "Fagus", "grandifolia",
  "R_902", 11.3, 10.6, 19.1,  "Fagus", "grandifolia"
)

results_many_many <- recruitMissingMatch(
  miss = missing_census_many_many, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_many_many, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 3, sd_dist = 1.5, inc_diam = 1.2, sd_inc_diam = 2,
  max_neg_diam = 2, max_rec_diam = 5, min_trust = 0.02)

results_many_many

Outcome: The stems form a “Many-to-Many” topology. The pairs dataframe provides trust scores for all combinations to facilitate final assignment.

Example 5: Shrinkage Penalty (max_neg_diam)

This tests the max_neg_diam penalty. A stem with an exact spatial and taxonomic match that exceeds the maximum allowable shrinkage is penalized.

missing_census_shrink <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 10.0, 10.0, 40.0,  "Pinus", "strobus"
)

recruit_census_shrink <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "R_901", 10.5, 10.2, 41.5,  "Pinus", "strobus",
  "R_902", 10.0, 10.0, 25.0,  "Pinus", "strobus"
)

results_shrink <- recruitMissingMatch(
  miss = missing_census_shrink, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_shrink, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5.0, sd_dist = 1.5, inc_diam = 1.0, sd_inc_diam = 2.0, 
  max_neg_diam = 2.0, max_rec_diam = 5.0, min_trust = 0.02
)

results_shrink

Outcome: R_901 is matched with T_101. R_902 receives a trust score below the 0.02 threshold due to the 0.01 shrinkage penalty multiplier and is routed to the true_recruits output.

Example 6: Spatial Cutoff (max_dist)

This tests the absolute spatial boundary defined by max_dist.

missing_census_dist <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 0.0,  0.0,  20.0,  "Acer",  "rubrum"
)

recruit_census_dist <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "R_901", 0.0,  5.1,  21.0,  "Acer",  "rubrum"
)

results_dist <- recruitMissingMatch(
  miss = missing_census_dist, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_dist, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5.0, sd_dist = 1.5, inc_diam = 1.0, sd_inc_diam = 2.0, 
  max_neg_diam = 2.0, max_rec_diam = 5.0, min_trust = 0.02
)

results_dist

Outcome: No matches are generated. R_901 is 5.1 m away, exceeding the 5.0 m max_dist limit, and is returned in the true_recruits data frame.

Example 7: Taxonomic Penalty

The algorithm applies a tiered multiplier for taxonomic mismatches: 1.0 for a complete match, 0.7 for a genus match, and 0.1 for a complete mismatch.

missing_census_taxa <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 15.0, 15.0, 30.0,  "Fagus", "grandifolia"
)

recruit_census_taxa <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,    ~species,
  "R_901", 15.0, 15.0, 31.0,  "Quercus", "rubra"
)

results_taxa <- recruitMissingMatch(
  miss = missing_census_taxa, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_taxa, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5.0, sd_dist = 1.5, inc_diam = 1.0, sd_inc_diam = 2.0, 
  max_neg_diam = 2.0, max_rec_diam = 5.0, min_trust = 0.15 
)

results_taxa

Outcome: The genus mismatch applies a 0.1 multiplier. With min_trust set to 0.15, the final trust score falls below the threshold, preventing a match.

Example 8: Disabling the Diameter Filter (max_rec_diam)

Setting max_rec_diam = 0.0 disables the initial size filter, forcing the evaluation of all stems regardless of diameter.

missing_census_tiny <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "T_101", 5.0,  5.0,  4.0,   "Acer",  "rubrum"
)

recruit_census_tiny <- tribble(
  ~id,     ~x,   ~y,   ~diam, ~genus,  ~species,
  "R_901", 5.2,  5.1,  4.5,   "Acer",  "rubrum"
)

# Execution 1: Default threshold
results_tiny_filtered <- recruitMissingMatch(
  miss = missing_census_tiny, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_tiny, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5.0, sd_dist = 1.5, inc_diam = 1.0, sd_inc_diam = 2.0, 
  max_neg_diam = 2.0, max_rec_diam = 5.0, min_trust = 0.02
)

results_tiny_filtered

# Execution 2: Filter disabled
results_tiny_forced <- recruitMissingMatch(
  miss = missing_census_tiny, 
  miss_id = "id", miss_x = "x", miss_y = "y", miss_diam = "diam", miss_genus = "genus", miss_species = "species",
  rec = recruit_census_tiny, 
  rec_id = "id", rec_x = "x", rec_y = "y", rec_diam = "diam", rec_genus = "genus", rec_species = "species",
  max_dist = 5.0, sd_dist = 1.5, inc_diam = 1.0, sd_inc_diam = 2.0, 
  max_neg_diam = 2.0, max_rec_diam = 0.0, min_trust = 0.02
)

results_tiny_forced

Outcome: Under default settings, R_901 is filtered as a true recruit. When max_rec_diam = 0.0, R_901 is evaluated and successfully matched to T_101.