1 Creating and printing a pig Simple functions used with a pig - - PDF document

1
SMART_READER_LITE
LIVE PREVIEW

1 Creating and printing a pig Simple functions used with a pig - - PDF document

System description A batch of slaughter pigs is simulated from insertion of the batch until The SimBatch simulation model slaughter. The pigs are fed ad libitum. Advanced Herd Management Sent to a Danish slaughter house. Anders Ringgaard


slide-1
SLIDE 1

1

Slide 1

The SimBatch simulation model

Advanced Herd Management Anders Ringgaard Kristensen

Slide 2

System description

A batch of slaughter pigs is simulated from insertion of the batch until slaughter. The pigs are fed ad libitum. Sent to a Danish slaughter house. The highest price per kg is obtained in the “optimal” interval from 70 to 84 kg slaughter weight. Delivery is decided on the basis of

  • bserved live weight.

The slaughter policy is decided by the farmer

Slide 3

Model type SimBatch is a

  • Dynamic (i.e. not static) time stepping model

with daily updating of the states of the pigs and the batch.

  • Mechanistic (i.e. not empirical) model since it

models a batch by modeling each individual pig.

  • Initially it is a deterministic model, but in the

exercise it is step by step transformed into a stochastic model.

  • Monte Carlo simulation model, because the

simulation mechanism is based on drawing random numbers from relevant distributions.

Slide 4

State of nature, Φ0

# State of nature information stateOfNature = list(initialAverage = 30, # Average weight at insertion initialStdDev = 2, # Std. deviation of weight at insertion averageDG = 0.900, # Average daily gain stdDevDG = 0.100, # Std. dev. of daily gain between pigs autocorrelationDG = 0.99, # Autocor. for individual daily gain diseaseRisk = 0.01, # Risk of disease any day diseaseEffectAve = 0.250, # Average effect of disease diseaseEffectStd = 0.030, # Standard deviation of disease effect diseaseDurationAv = 5, # Average duration of disease liveToSlaughterAv = 0.76, # Live w. to slaughter w. conversion liveToSlaughterStd = 1.4, # Standard deviation of conversion lowerOptimal = 70, # Lower bound of opt. slaughter weight upperOptimal = 84, # Upper bound of opt. slaughter weight deliverEvery = 7) # Days btw. del.(7 = weekly, 1 = daily) Slide 5

Decision strategy, Θ

# Decision strategy information decisionStrategy = list(thresholdWeight = 100, # Intended live weight at slaughter minimumDelivery = 10, # Minimum number of pigs to send to maximumUnderWeight = 2, # Acc this number of pigs bel. Thresh. minimumStock = 50, # Minimum acc. stock: Send all if below maximumAge = 90) # Maximum age before slaughter # Strategy with completely individual delivery individualStrategy = list(thresholdWeight = 100, # Intended live weight at slaughter minimumDelivery = 1, # Minimum number of pigs to send to maximumUnderWeight = 0, # Acc this number of pigs bel. thres minimumStock = 0, # Minimum acc. stock: Send all if bel maximumAge = 9999) # Maximum age before slaughter # Strategy with slaughter at a given age ageStrategy = list(thresholdWeight = 9999, # Intended live weight at slaughter minimumDelivery = 9999, # Minimum number of pigs to send to maximumUnderWeight = 0, # Acc this number of pigs bel. thres minimumStock = 0, # Minimum acc. stock: Send all if bel maximumAge = 90) # Maximum age before slaughter Slide 6

Modeling a pig

A pig is defined by 4 state variables

  • Live weight in kg
  • The gain today in kg (if healthy)
  • An integer, i, defining the status of the pig:
  • i = -2: The pig has already been slaughtered
  • i = -1: The pig is healthy and still present
  • i ≥ 0: The pig is diseased and will stay diseased i

days yet

  • The effect of disease (if any) on daily gain

In SimBatch, a pig is represented on a given day by a vector with 4 elements. The state of the pig is updated every day.

slide-2
SLIDE 2

2

Slide 7

Creating and printing a pig A valid pig object: > pig = c(45.0, 0.956, 3, 0.123) > pig [ 1] 45.000 0.956 3.000 0.123 A pig with a live weight of 45 kg, a gain until tomorrow of 956 g (if healthy). It is, however diseased and will stay diseased for 3 days. The effect of disease is a reduction in daily gain of 123 g

Slide 8

Simple functions used with a pig object

getStartW eight( stateOfNature) : A start weight (in kg) is returned. Initially, it only returns the initial average from the state of nature. getI nitialDG( stateOfNature) : An initial daily gain (in kg) for a pig is

  • returned. Initially, it only returns the average daily gain from the state of

nature. isPigDiseased( stateOfNature) : Returns a logical value. If the pig gets diseased it returns TRUE, otherwise FALSE. Initially, it always returns FALSE. getDiseaseEffect( stateOfNature) : Returns the reduction (in kg) of daily gain during a disease period. Initially it only returns the average disease effect from the state of nature. getDiseaseDuration( stateOfNature) : Returns the duration (in days) of a new disease. Initially, it only returns the average duration from the state of nature. getSlaughterW eight( liveW eight, stateOfNature) : Transforms the specified input parameter liveWeight to slaughter weight (in kg). Initially it just multiplies the live weight with the conversion ratio from the state of nature. getObservedW eight( liveW eight) : Returns the observed live weight of a pig with the specified true live weight. Initially, it just returns the true live weight (ignoring the observation error). updateDG( oldDG, stateOfNature) : Returns a new daily gain for the next 24 hours from the specified “old” daily gain (from the previous 24 hours). Initially, it just returns the old value (without any transform).

Slide 9

Function creating a piglet (don’t change)

# Create a "30 kg" piglet at random createPiglet = function(son) { piglet = rep(NA, 4) piglet[ 1] = getStartWeight(son) # Draw initial live weight piglet[ 2] = getInitialDG(son) # Draw initial daily gain piglet[ 3] = -1 # Number of days remaining if diseased # If "-1": Healthy. If "-2": Slaughtered piglet[ 4] = 0 # Disease effect on daily gain return(piglet) }

Slide 10

Function updating a pig (don’t change)

# Updates the state of a pig updatePig = function(pig, son) { pigNew = rep(NA, 4) If (pig[3] > -2) { # Is the pig still present? pigNew[1] = pig[1] + pig[2] - pig[4] # Yes, calculate new weight pigNew[2] = updateDG(pig[2], son) # Update daily gain if (pig[3] == -1) { # Is the pig healthy now? if (isPigDiseased(son)) { # Yes it was. Should it change state? pigNew[3] = getDiseaseDuration(son) # Yes it should. Draw duration of disease pigNew[4] = getDiseaseEffect(son) # Draw effect of disease } else { # No, the pig remains healthy pigNew[3] = -1 # Define it as healthy pigNew[4] = 0 # No effect of disease } } else { # No, the pig is already diseased if (pig[3] > 0) { # Has it still at least one day left of disease? pigNew[3] = pig[3] - 1; # Yes, reduce the number of days left by one pigNew[4] = pig[4] # Keep the disease effect on daily gain } else { # No, it is the last day of the disease pigNew[3] = -1; # Change state to healthy pigNew[4] = 0 # No disease effect next time } } } else { # The pig is already slaughtered pigNew = pig # Leave it as it is } return(pigNew) } Slide 11

Representation of a batch on a given day

Just an array of pigs. Exam ple of a batch w ith 5 piglets:

> pig1 = createPiglet(stateOfNature) > pig2 = createPiglet(stateOfNature) > pig3 = createPiglet(stateOfNature) > pig4 = createPiglet(stateOfNature) > pig5 = createPiglet(stateOfNature) > batchDay = array(NA, dim = c(5, length(pig1))) > batchDay[ 1,] = pig1 > batchDay[ 2,] = pig2 > batchDay[ 3,] = pig3 > batchDay[ 4,] = pig4 > batchDay[ 5,] = pig5 > batchDay [ ,1] [ ,2] [ ,3] [ ,4] [ 1,] 29.77745 0.8219607 -1 0 [ 2,] 27.78206 0.9373041 -1 0 [ 3,] 27.52593 0.8217774 -1 0 [ 4,] 31.96445 0.9472012 -1 0 [ 5,] 30.67659 0.8425552 -1 0

Slide 12

Function updating a batch (don’t change)

# Updates the state of a batch updateBatch = function(batchDay, son) { newBatchDay = array(NA, dim = c(length(batchDay[ ,1] ), length(batchDay[ 1,] ))) for (i in 1: length(batchDay[ ,1] )) { # Iterate over pigs pigDay = batchDay[ i,] # Take the i’th pig pigDay = updatePig(pigDay, son) # Update it newBatchDay[ i,] = pigDay # Insert it in the updated batch } return(newBatchDay) # Return the updated batch }

Function counting pigs left in a batch (don’t change)

# Counts the number of pigs still present in the batch (not slaughtered) pigsLeft = function(batchDay) { left = 0 for (i in 1: length(batchDay[ ,1] )) { if (batchDay[ i,3] > -2) { left = left + 1 } } return(left) }
slide-3
SLIDE 3

3

Slide 13

Applying decision strategy to batch on a given day

# Applies the decision strategy on the batch a given day selectForSlaughter = function(batchDay, age, dec) { correctedDay = array(NA, dim= c(length(batchDay[ ,1] ), length(batchDay[ 1,] )))
  • bservedDay = rep(NA, length(batchDay[ ,1] ))
for (i in 1: length(batchDay[ ,1] )) { correctedDay[ i,] = batchDay[ i,] # Copy pigs to new array
  • bservedDay[ i] = getObservedWeight(batchDay[ i,1] )
# Get observed weights } leftPigs = pigsLeft(batchDay) # How many pigs are still present if (leftPigs < dec$minimumStock | | # Fewer than minimum? age > dec$maximumAge) { # Maximum age exceeded? for (i in 1: length(batchDay[ ,1] )) { # Yes, correctedDay[ i,3] = -2 # ... cull the rest } } readyI ndexes = which(observedDay > = dec$thresholdWeight) # Find the pigs with observed weights exceeding threshold readyCount = length(readyI ndexes) # How many were there if (readyCount < dec$minimumDelivery) { # Are there two few? if (readyCount > dec$minimumDelivery - dec$maximumUnderWeight) { # Yes, but are there so many that we slaughter anyway? sortedVector = sort(observedDay, decreasing = TRUE) # Yes, find the biggest pigs realThres = sortedVector[ dec$minimumDelivery] # Find the real threshold readyIndexes = which(observedDay > = realThres) # Pick all pigs exceeding the real threshold for (i in readyI ndexes) { correctedDay[ i,3] = -2 # Mark them as slaughtered } } } else { # No, there were enough pigs for (i in readyIndexes) { # Pick them correctedDay[ i,3] = -2 # Mark them as slaughtered } } return(correctedDay) # Return the pigs with the slaughtered marked } Slide 14

A simulation run – function simulateBatch (don’t change)

# Simulates a batch with the specified number of pigs over the specified number of days simulateBatch = function(pigs, days, son, dec) { piglet = createPiglet(son) # Create a piglet (just to know its dimension ...) batchStructure = array(NA, dim= c(pigs, days, length(piglet))) # Create a data structure for the results for (i in 1: pigs) { # Create "pigs" piglets piglet = createPiglet(son) # The i’th piglet batchStructure[ i,1,] = piglet # Insert it } for (n in 2: days) { # Iterate over the simulation period day by day batchStructure[ ,n,] = updateBatch(batchStructure[ ,n-1,] , son) # if (n % % son$deliverEvery = = 0) { # Is it a delivery day? batchStructure[ ,n,] = selectForSlaughter(batchStructure[ ,n,] , n, dec) # Yes, apply decision strategy } } return(batchStructure) # Return full data structure }

Slide 15

Plotting weight curves from a simulation run

# Plot growth curves of pigs plotBatchWeights = function(batch, skip) { pigWeights = batch[ 1,,] days = 1: length(pigWeights[ ,1] ) plot(days, pigWeights[ ,1] , type= "l", ylim= c(25, 110), col= "white", xlab= "Days", ylab= "Kg live weight") for (i in 1: length(batch[ ,1,1] )) { pigLife = batch[ i,,] liveIndexes = which(pigLife[ ,3] > -2) pigWeights = pigLife[ liveIndexes,1] if ((i-1) % % skip = = 0) { lines(pigWeights) } } }

Slide 16

Plottet weight curves from a simulation run

Slide 17

Functions returning output variables

batch = simulateBatch(100, 120, stateOfNature, decisionStrategy)

  • Runs a simulation with 100 pigs over 120 days using the specified state of

nature and the specified decision strategy.

plotBatchWeights(batch, 1)

  • Plots the resulting weight curves

deliveredInOptimalInterval(batch, stateOfNature)

  • Calculate (and print) the number of pigs delivered in the optimal interval

dailyGain(batch)

  • Calculates the true average daily gain for the batch