1
votes

I have a simulated data set where each row is an individual bird, and I've written some functions that determine whether each individual lives or dies. The conditions for whether or not an individual survives is based off of its age class (AHY or HY) and its sex (M or F). I created a function for each age/sex combination, and use pmap_chr inside of mutate/case_when, which should fill in a column called status. In my code this provides a value of 'live' or 'die.' Here's an abbreviated version of my dataset:

library(tidyverse)
agents <- structure(list(id = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 
13, 14, 15, 16, 17, 18), sex = c("F", "F", "F", "F", "F", "M", 
"M", "M", "M", "M", "F", "F", "M", "M", "M", "M", "M", "M"), 
    class = c("AHY", "AHY", "AHY", "AHY", "AHY", "AHY", "AHY", 
    "AHY", "AHY", "AHY", "HY", "HY", "HY", "HY", "HY", "HY", 
    "HY", "HY"), hDateCtr = c(-0.84852029241304, 0.558881154137435, 
    -0.909711659654365, 1.21158907137824, -0.56296057862019, 
    -0.0938267631033649, -1.54202245448139, -0.216209497586015, 
    1.33397180586089, 1.06880921448181, -0.935414346693485, -0.935414346693485, 
    -0.935414346693485, -0.935414346693485, 0.935414346693485, 
    0.935414346693485, 0.935414346693485, 0.935414346693485), 
    aDateCtr = c(-1.13245629117638, 1.13245629117638, -0.490731059509763, 
    1.13245629117638, -0.641725231666613, 1.13245629117638, -1.13245629117638, 
    1.13245629117638, -0.490731059509763, -0.641725231666613, 
    NA, NA, NA, NA, NA, NA, NA, NA), selfOrig = c("imm", "imm", 
    "imm", "imm", "imm", "imm", "imm", "imm", "imm", "imm", "local", 
    "local", "local", "local", "local", "local", "local", "local"
    ), sameSexOrig = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
    "imm", "imm", "imm", "imm", "imm", "imm", "imm", "imm"), 
    success = c(TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, 
    TRUE, FALSE, FALSE, NA, NA, NA, NA, NA, NA, NA, NA), paired = c(TRUE, 
    TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, 
    FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE)), class = c("tbl_df", 
"tbl", "data.frame"), row.names = c(NA, -18L))

# A tibble: 18 x 9
      id sex   class hDateCtr aDateCtr selfOrig sameSexOrig success paired
   <dbl> <chr> <chr>    <dbl>    <dbl> <chr>    <chr>       <lgl>   <lgl> 
 1     1 F     AHY    -0.849    -1.13  imm      NA          TRUE    TRUE  
 2     2 F     AHY     0.559     1.13  imm      NA          TRUE    TRUE  
 3     3 F     AHY    -0.910    -0.491 imm      NA          FALSE   TRUE  
 4     4 F     AHY     1.21      1.13  imm      NA          FALSE   TRUE  
 5     5 F     AHY    -0.563    -0.642 imm      NA          FALSE   TRUE  
 6     6 M     AHY    -0.0938    1.13  imm      NA          FALSE   TRUE  
 7     7 M     AHY    -1.54     -1.13  imm      NA          TRUE    TRUE  
 8     8 M     AHY    -0.216     1.13  imm      NA          TRUE    TRUE  
 9     9 M     AHY     1.33     -0.491 imm      NA          FALSE   TRUE  
10    10 M     AHY     1.07     -0.642 imm      NA          FALSE   TRUE  
11    11 F     HY     -0.935    NA     local    imm         NA      FALSE 
12    12 F     HY     -0.935    NA     local    imm         NA      FALSE 
13    13 M     HY     -0.935    NA     local    imm         NA      FALSE 
14    14 M     HY     -0.935    NA     local    imm         NA      FALSE 
15    15 M     HY      0.935    NA     local    imm         NA      FALSE 
16    16 M     HY      0.935    NA     local    imm         NA      FALSE 
17    17 M     HY      0.935    NA     local    imm         NA      FALSE 
18    18 M     HY      0.935    NA     local    imm         NA      FALSE

Here's an example of the mortality functions I've written that go into pmap_chr. These all work fine if I run the code below on a data set that has a single age class or sex:

hDateEffect <- TRUE
winterTemp <- -3

# hatchling mortality -----------------------------------------------------

hatchMortInt <- -4.67
hatchMortIntSD <- 0.39
hatchMortBeta1 <- 0.6
hatchMortBeta1SD <- 0.27
hatchMortBeta2 <- 1.12
hatchMortBeta2SD <- 0.36
hatchMortBeta3 <- -0.3
hatchMortBeta3SD <- 0.16
hatchMortBeta4 <- -0.3
hatchMortBeta4SD <- 0.16


# male mortality ----------------------------------------------------------

maleMortInt <- -2.09
maleMortIntSD <- 0.32
maleMortBeta1 <- 0.81
maleMortBeta1SD <- 0.34
maleMortBeta2 <- -1.36
maleMortBeta2SD <- 0.84
maleMortBeta3 <- 1.67
maleMortBeta3SD <- 0.32


# female mortality --------------------------------------------------------

femMortInt <- -0.93
femMortIntSD <- 0.87
femMortBeta1 <- 1.59
femMortBeta1SD <- 0.35
femMortBeta2 <- -1.77
femMortBeta2SD <- 0.78

# hatch-year female 
HY_female_mortality <- function(hDateCtr, sameSexOrig, ...) {
  intercept <- rnorm(1, hatchMortInt, hatchMortIntSD)
  beta2 <- rnorm(1, hatchMortBeta2, hatchMortBeta2SD)
  beta4 <- rnorm(1, hatchMortBeta4, hatchMortBeta4SD)
  if(hDateEffect == TRUE) {
    beta3 <- rnorm(1, hatchMortBeta3, hatchMortBeta3SD)
  } else {
    beta3 <- 0
  }
  if (sameSexOrig == 'local') {
    linSurv <- intercept + beta2 + (beta3 * hDateCtr) + (beta4 * winterTemp)
  } else {
    linSurv <- intercept + (beta3 * hDateCtr) + (beta4 * winterTemp)
  }
  probSurv <- plogis(linSurv)
  randDraw <- runif(1, 0, 1)
  if (randDraw > probSurv) {
    val <- 'die'
    return(val)
  } else {
    val <- 'live'
    return(val)
  }
}

# hatch-year male
HY_male_mortality <- function(hDateCtr, sameSexOrig, ...) {
  intercept <- rnorm(1, hatchMortInt, hatchMortIntSD)
  beta1 <- rnorm(1, hatchMortBeta1, hatchMortBeta1SD) 
  beta2 <- rnorm(1, hatchMortBeta2, hatchMortBeta2SD)
  beta4 <- rnorm(1, hatchMortBeta4, hatchMortBeta4SD)
  if(hDateEffect == TRUE) {
    beta3 <- rnorm(1, hatchMortBeta3, hatchMortBeta3SD)
  } else {
    beta3 <- 0
  }
  if (sameSexOrig == 'local') {
    linSurv <- intercept + beta1 + beta2 + (beta3 * hDateCtr) + (beta4 * winterTemp)
  } else {
    linSurv <- intercept + beta1 + (beta3 * hDateCtr) + (beta4 * winterTemp)
  }
  probSurv <- plogis(linSurv)
  randDraw <- runif(1, 0, 1)
  if (randDraw > probSurv) {
    val <- 'die'
    return(val)
  } else {
    val <- 'live'
    return(val)
  }
}

# after-hatch-year mortality functions 
# after-hatch-year male 
AHY_male_mortality <- function(aDateCtr, success, selfOrig, ...) {
  intercept <- rnorm(1, maleMortInt, maleMortIntSD)
  beta1 <- rnorm(1, maleMortBeta1, maleMortBeta1SD) 
  beta3 <- rnorm(1, maleMortBeta3, maleMortBeta3SD)
  if(hDateEffect == TRUE) {
    beta2 <- rnorm(1, hatchMortBeta3, hatchMortBeta3SD)
  } else {
    beta2 <- 0
  }
  if (success == TRUE) {
    linSurv <- intercept + beta1 + (beta2 * aDateCtr)
  } else {
    linSurv <- intercept + (beta2 * aDateCtr)
  }
  if (selfOrig == 'local') {
    linSurv <- linSurv + beta3
  } else {
    linSurv <- linSurv
  }
  probSurv <- plogis(linSurv)
  randDraw <- runif(1, 0, 1)
  if (randDraw > probSurv) {
    val <- 'die'
    return(val)
  } else {
    val <- 'live'
    return(val)
  }
}

# after-hatch-year female 
AHY_female_mortality <- function(aDateCtr, success, ...) {
  intercept <- rnorm(1, femMortInt, femMortIntSD)
  beta1 <- rnorm(1, femMortBeta1, femMortBeta1SD) 
  beta2 <- rnorm(1, femMortBeta2, femMortBeta2SD)
  if (success == TRUE) {
    linSurv <- intercept + beta1 + (beta2 * aDateCtr)
  } else {
    linSurv <- intercept + (beta2 * aDateCtr)
  }
  probSurv <- plogis(linSurv)
  randDraw <- runif(1, 0, 1)
  if (randDraw > probSurv) {
    val <- 'die'
  } else {
    val <- 'live'
  }
  return(val)
}

Here's the pmap_chr part, which doesn't work across all of the age and sex class combinations:

agents %>%
  mutate(
    status = case_when(
      class == 'HY' & sex == 'F' ~ pmap_chr(., HY_female_mortality),
      class == 'HY' & sex == 'M' ~ pmap_chr(., HY_male_mortality),
      class == 'AHY' & sex == 'M' ~ pmap_chr(., AHY_male_mortality),
      class == 'AHY' & sex == 'F' ~ pmap_chr(., AHY_female_mortality)
    )
  )

But if I do this same thing instead for another logical called 'success' (so if (success == TRUE)) which is the one I actually need the condition to be based, it produces an error:

Error in mutate_impl(.data, dots) : 
  Evaluation error: missing value where TRUE/FALSE needed.

I'm at a loss for why these functions work separately, but not on the entire dataset containing all age and sex classes. I have examples of different processes (reproduction, immigration) where I do a similar thing (take the data set, write functions that are used inside pmap, which in turn are inside case_when and mutate).

2
pmap_chr(agents, HY_female_mortality) throws an error for me, which is ultimately what you're calling in your case_when statement. I think it's due to an NA value being to your sameSexOrig value and the design of your user-defined functions.zack
@zack it seems to work though if there are just HY females in the data set though (this works): agents %>% filter(class == 'HY' & sex == 'F') %>% mutate( status = case_when( class == 'HY' & sex == 'F' ~ pmap_chr(., HY_female_mortality) ) )Jason
yes, that is because sameSexOrig does not contain NA values for those observations.zack
@zack ok see what you mean. any idea why the the AHY functions, which have NAs for sameSexOrig but don't rely on that variable won't work?Jason

2 Answers

0
votes

Splitting by age class and taking a similar approach solves the problem, still not sure why it won't work the original way though...

agents <- agents %>%
    split(.$class) 

  agents$HY <- agents$HY %>%
    mutate(
      status = case_when(
        sex == 'F' ~ pmap_chr(., HY_female_mortality),
        sex == 'M' ~ pmap_chr(., HY_male_mortality)
      )
    )

  agents$AHY <- agents$AHY %>%
    mutate(
      status = case_when(
        sex == 'F' ~ pmap_chr(., AHY_female_mortality),
        sex == 'M' ~ pmap_chr(., AHY_male_mortality)
      )
    )

  agents <- agents %>%
    bind_rows()
0
votes

In response to A. Suliman's comment, I changed the functions so you can see that they give a particular value for each age class and sex:

# hatch-year female 

HY_female_mortality <- function(hDateCtr, sameSexOrig, ...) {
  if(hDateEffect == TRUE) {
    val <- 'hatch effect on'
  } else {
    val <- 'hatch effect off'
  }
  if (sameSexOrig == 'local') {
    val <- paste0(val, ' and local')
  } else {
    val <- paste0(val, ' and immigrant')
  }
    return(paste0(val, ' and female HY'))
}

# hatch-year male
HY_male_mortality <- function(hDateCtr, sameSexOrig, ...) {
  if(hDateEffect == TRUE) {
    val <- 'hatch effect on'
  } else {
    val <- 'hatch effect off'
  }
  if (sameSexOrig == 'local') {
    val <- paste0(val, ' and local')
  } else {
    val <- paste0(val, ' and immigrant')
  }
  return(paste0(val, ' and male HY'))
}

# after-hatch-year mortality functions 
# after-hatch-year male 
AHY_male_mortality <- function(aDateCtr, success, selfOrig, ...) {
  if(hDateEffect == TRUE) {
    val <- 'hatch effect on'
  } else {
    val <- 'hatch effect off'
  }
  if (success == TRUE) {
    val <- paste0(val, ' and successful')
  } else {
    val <- paste0(val, ' and failed')
  }
  if (selfOrig == 'local') {
    val <- paste0(val, ' and local')
  } else {
    val <- paste0(val, ' and immigrant')
  }
  return(paste0(val, ' and male AHY'))
}

# after-hatch-year female 
AHY_female_mortality <- function(aDateCtr, success, ...) {
  if(hDateEffect == TRUE) {
    val <- 'hatch effect on'
  } else {
    val <- 'hatch effect off'
  }
  if (success == TRUE) {
    val <- paste0(val, ' and successful')
  } else {
    val <- paste0(val, ' and failed')
  }
  return(paste0(val, ' and female AHY'))
}

agents <- agents %>%
  split(.$class) 

agents$HY %>%
  mutate(
    status = case_when(
      sex == 'F' ~ pmap_chr(., HY_female_mortality),
      sex == 'M' ~ pmap_chr(., HY_male_mortality)
    )
  )

agents$AHY %>%
  mutate(
    status = case_when(
      sex == 'F' ~ pmap_chr(., AHY_female_mortality),
      sex == 'M' ~ pmap_chr(., AHY_male_mortality)
    )
  )

Is this not behaving correctly?

> agents$HY
# A tibble: 8 x 10
     id sex   class hDateCtr aDateCtr selfOrig sameSexOrig success paired status                                     
  <dbl> <chr> <chr>    <dbl>    <dbl> <chr>    <chr>       <lgl>   <lgl>  <chr>                                      
1    11 F     HY      -0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and female HY
2    12 F     HY      -0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and female HY
3    13 M     HY      -0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY  
4    14 M     HY      -0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY  
5    15 M     HY       0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY  
6    16 M     HY       0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY  
7    17 M     HY       0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY  
8    18 M     HY       0.935       NA local    imm         NA      FALSE  hatch effect on and immigrant and male HY

> agents$AHY
# A tibble: 10 x 10
      id sex   class hDateCtr aDateCtr selfOrig sameSexOrig success paired status                                                   
   <dbl> <chr> <chr>    <dbl>    <dbl> <chr>    <chr>       <lgl>   <lgl>  <chr>                                                    
 1     1 F     AHY    -0.849    -1.13  imm      NA          TRUE    TRUE   hatch effect on and successful and female AHY            
 2     2 F     AHY     0.559     1.13  imm      NA          TRUE    TRUE   hatch effect on and successful and female AHY            
 3     3 F     AHY    -0.910    -0.491 imm      NA          FALSE   TRUE   hatch effect on and failed and female AHY                
 4     4 F     AHY     1.21      1.13  imm      NA          FALSE   TRUE   hatch effect on and failed and female AHY                
 5     5 F     AHY    -0.563    -0.642 imm      NA          FALSE   TRUE   hatch effect on and failed and female AHY                
 6     6 M     AHY    -0.0938    1.13  imm      NA          FALSE   TRUE   hatch effect on and failed and immigrant and male AHY    
 7     7 M     AHY    -1.54     -1.13  imm      NA          TRUE    TRUE   hatch effect on and successful and immigrant and male AHY
 8     8 M     AHY    -0.216     1.13  imm      NA          TRUE    TRUE   hatch effect on and successful and immigrant and male AHY
 9     9 M     AHY     1.33     -0.491 imm      NA          FALSE   TRUE   hatch effect on and failed and immigrant and male AHY    
10    10 M     AHY     1.07     -0.642 imm      NA          FALSE   TRUE   hatch effect on and failed and immigrant and male AHY