RBFY Pen-Level DMI Spatial Analysis

Description

This report characterizes spatial patterns of dry matter intake (DMI) across pens at Royal Beef Feedyard using operational data from June 2022 to current. The analytical unit is the pen-week. To be included, a pen-week was required to represent a designated feeding pen with a weekly average of at least 30 animals, have at least five days of feeding records during the week, and maintain consistent lot occupancy throughout the week. Hospital, receiving, processing, sort, and shipping pens were excluded.

Descriptive DMI heat map

June 2022 to current

MBW = metabolic body weight (BW^0.75)

summaries   = transpose(pen_summaries_ojs)
coords      = transpose(pen_coords_ojs)
metric_info = transpose(metric_meta_ojs)

// Dropdowns — metric first, defaults to first level (DMI, lb/hd/day)
viewof metric_sel = Inputs.select(metric_levels_ojs, {label: "Metric"})
viewof sex_sel    = Inputs.select(sex_levels_ojs,    {label: "Sex"})
viewof weight_sel = Inputs.select(weight_levels_ojs, {label: "Weight"})
viewof season_sel = Inputs.select(season_levels_ojs, {label: "Season"})

selected = summaries.filter(d =>
  d.metric === metric_sel &&
  d.sex    === sex_sel    &&
  d.weight === weight_sel &&
  d.season === season_sel
)

L = require("leaflet@1.9.4/dist/leaflet-src.js")
html`<link rel="stylesheet" href="https://unpkg.com/leaflet@1.9.4/dist/leaflet.css">`

{
  const geojson = JSON.parse(pen_geojson_ojs)

  function makeLookup(sel) {
    return Object.fromEntries(sel.map(d => [d.PEN, d]))
  }

  function getColor(val, minV, maxV) {
    if (val === null || val === undefined || isNaN(val)) return "#808080"
    const t = Math.max(0, Math.min(1, (val - minV) / (maxV - minV)))
    const stops = [
      [215,  48,  39],
      [252, 141,  89],
      [254, 224, 139],
      [255, 255, 191],
      [217, 239, 139],
      [145, 207,  96],
      [ 26, 152,  80]
    ]
    const scaled = t * (stops.length - 1)
    const lo     = Math.floor(scaled)
    const hi     = Math.min(lo + 1, stops.length - 1)
    const frac   = scaled - lo
    const r = Math.round(stops[lo][0] + frac * (stops[hi][0] - stops[lo][0]))
    const g = Math.round(stops[lo][1] + frac * (stops[hi][1] - stops[lo][1]))
    const b = Math.round(stops[lo][2] + frac * (stops[hi][2] - stops[lo][2]))
    return `rgb(${r},${g},${b})`
  }

function quantileSorted(values, p) {
    if (values.length === 0) return null
    const sorted = [...values].sort((a, b) => a - b)
    const idx    = (sorted.length - 1) * p
    const lo     = Math.floor(idx)
    const hi     = Math.ceil(idx)
    if (lo === hi) return sorted[lo]
    return sorted[lo] + (idx - lo) * (sorted[hi] - sorted[lo])
  }

  const div = document.createElement("div")
  div.style.height = "600px"
  div.style.width  = "100%"
  yield div

  const lats   = coords.map(d => d.lat)
  const lons   = coords.map(d => d.lon)
  const minLat = Math.min(...lats), maxLat = Math.max(...lats)
  const minLon = Math.min(...lons), maxLon = Math.max(...lons)

  const map = L.map(div, { minZoom: 15.25 })
    .fitBounds([[minLat - 0.01, minLon - 0.01],
                [maxLat + 0.01, maxLon + 0.01]])

  map.setMaxBounds([
    [minLat - 0.003, minLon - 0.003],
    [maxLat + 0.003, maxLon + 0.003]
  ])

  L.tileLayer(
    "https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{z}/{y}/{x}",
    { attribution: "Tiles © Esri" }
  ).addTo(map)

  let polyLayer = L.geoJSON(null).addTo(map)

  function updateMap(sel) {
    // Metadata for the selected metric (units label + decimals)
    const meta     = metric_info.find(m => m.label === metric_sel)
    const decimals = meta ? meta.decimals : 2
    // Units string for legend/popups: strip "DMI, " prefix from the label
    const units    = meta ? meta.label.replace(/^DMI,\s*/, "") : ""

    const lkp  = makeLookup(sel)
    const vals = Object.values(lkp)
      .map(d => d.mean_val)
      .filter(v => v !== null && v !== undefined && !isNaN(v))
    // Percentile-capped scale bounds (1st/99th) so outlier pens don't
    // compress the color range; values beyond the caps take the end colors.
    const minV = vals.length > 0 ? quantileSorted(vals, 0.005) : 0
    const maxV = vals.length > 0 ? quantileSorted(vals, 0.995) : 1

    polyLayer.clearLayers()
    polyLayer = L.geoJSON(geojson, {
      style: feature => {
        const d   = lkp[feature.properties.pen_id]
        const val = (d && d.mean_val !== null && !isNaN(d.mean_val))
                      ? d.mean_val : null
        return {
          fillColor:   getColor(val, minV, maxV),
          fillOpacity: 0.6,
          color:       "white",
          weight:      1
        }
      },
      onEachFeature: (feature, layer) => {
        const penId   = feature.properties.pen_id
        const d       = lkp[penId]
        const hasData = d && d.mean_val !== null && !isNaN(d.mean_val)
        const fmt     = v => (v !== null && v !== undefined && !isNaN(v))
                               ? `${Number(v).toFixed(decimals)} ${units}` : "N/A"
        layer.bindTooltip(
          hasData ? `Pen ${penId}: ${fmt(d.mean_val)}`
                  : `Pen ${penId}: no data`,
          { direction: "auto" }
        )
        layer.bindPopup(
          hasData
            ? `<b>Pen ${penId}</b><br>
               Mean: ${fmt(d.mean_val)}<br>
               Median: ${fmt(d.median_val)}<br>
               SD: ${fmt(d.sd_val)}<br>
               Pen-weeks: ${d.n_pen_weeks}`
            : `<b>Pen ${penId}</b><br>No data for this selection`
        )
      }
    }).addTo(map)

    const legend = L.control({ position: "bottomright" })
    legend.onAdd = function() {
      const el = L.DomUtil.create("div")
      el.style.cssText = `
        background: white; padding: 8px 10px; border-radius: 4px;
        box-shadow: 0 1px 4px rgba(0,0,0,0.3); font-size: 14px; line-height: 1.6;
      `
      const steps  = 100
      const height = 150
      let bars = ""
      for (let i = steps; i >= 0; i--) {
        const t   = i / steps
        const val = minV + t * (maxV - minV)
        bars += `<div style="width:16px; height:${height/steps}px;
                             background:${getColor(val, minV, maxV)};"></div>`
      }
      const nTicks = 5
      let ticks = ""
      for (let i = nTicks; i >= 0; i--) {
        const t   = i / nTicks
        const val = minV + t * (maxV - minV)
        const top = ((nTicks - i) / nTicks) * height
        const tickLabel = (i === nTicks)      ? `≤ ${val.toFixed(decimals)}`
                        : (i === 0)           ? `≥ ${val.toFixed(decimals)}`
                        :                        val.toFixed(decimals)
        ticks += `<div style="position:absolute; top:${top}px; left:22px;
                              font-size:14px; white-space:nowrap;
                              transform:translateY(-50%);">
                    ${tickLabel}
                  </div>`
      }
      el.innerHTML = `
        <b style="display:block; margin-bottom:4px;">${metric_sel}</b>
        <div style="position:relative; display:flex; align-items:stretch;">
          <div style="display:flex; flex-direction:column;">${bars}</div>
          <div style="position:relative; width:64px;">${ticks}</div>
        </div>
      `
      return el
    }
    if (map._legend) map.removeControl(map._legend)
    map._legend = legend
    legend.addTo(map)
  }

  invalidation.then(() => map.remove())
  updateMap(selected)
}

The color scale anchors to the minimum and maximum values within the currently selected combination (i.e., the legend scale changes with each change to filterable parameters).

Baseline mixed model

A multivariable linear mixed model was fit to estimate each pen’s intake after adjusting for the major non-spatial sources of variation in DMI. The fixed effects below were included as adjustment terms so that differences in cattle type, ration, health, and season are not mistaken for pen-driven differences in intake:

• Sex, breed, background, risk classification, and origin
• Current estimated BW
• Days on feed
• Ration NEg
• Natural-program status
• Cumulative mortality rate
• Weekly morbidity rate
• Days since last implant
• Month of year

Lot and year were included as random effects. Lot accounts for intake being heavily shaped by the specific cohort of cattle occupying a pen, and year accounts for broad temporal shifts from year-to-year variability; including them keeps cohort and time-period effects from leaking into the pen estimates. Pen was also included as a random effect, and its estimated values are mapped - representing each pen’s tendency to feed above or below expectation once all of the above has been accounted for.

Model-adjusted heat map

The color scale shows each pen’s adjusted intake as a deviation from the yard average - these values are adjusted for cattle characteristics, management, and season, but not for physical pen attributes such as size, location, or bunk space. A pen that remains consistently above or below average therefore indicates a location effect worth further examination against physical pen features.

These adjusted values control for measured differences in cattle, ration, health, and season, but they cannot account for unmeasured cattle characteristics that influence intake. Where cattle are placed non-randomly - for example, if lighter cattle are concentrated in certain pens - any intake differences not captured by the measured variables (e.g., genetics of the lighter cattle) will be attributed to the pen rather than to the cattle. A consistently above- or below-average pen should therefore be read as a flag for further investigation rather than as a confirmed location effect.