Estimating Population Size using RMark

Question

I am estimating the population size of a species using the combination of "marked" and "RMark" analysis. I am beginner to R and I follow someone's script here. My dataset are from capture-recapture obs. and exactly looks like dipper data on "marked" package but with fewer capturing events.

There are 3 steps:

1. Finding the best model using "marked" and getting P, phi, pent and N parameters.
2. Estimating the population size using the parameters on step 1.
3. Estimating another population size with "RMark" package and the best model on step 1 (the script i followed said it will estimate with 95% confidence level)

Here I attached the full script to make everyone easier to help me:

Step 1

library(marked)
moth<-read.csv("moth.eps1.csv")
View(moth)
moth.js.proc <- process.data(moth, model = "JS", groups = "sex")
moth.js.ddl <- make.design.data(moth.js.proc)
fit.js.moth.models <- function(){
  Phi.dot <- list(formula=~1)
  Phi.time <- list(formula=~time)
  p.dot <- list(formula=~1)
  pent.time <- list(formula=~time)
  pent.sex <- list(formula=~sex)
  pent.dot <- list(formula=~1)
  N.sex <- list(formula=~sex)
  N.dot <- list(formula=~1)
  cml <- create.model.list(c("Phi","p", "pent", "N"))
  results <- crm.wrapper(cml, data = moth.js.proc, ddl = moth.js.ddl,
                         external = FALSE, accumulate = FALSE, hessian = TRUE)
  
  return(results)
}
moth.js.models <- fit.js.moth.models()
moth.js.models
moth.js.models[[1]]
moth.js.predicted <- predict(moth.js.models[[1]]) # [[1]] just calls the model row according to the model table.
moth.js.predicted #run better from the console

Step 2

N.derived <- data.frame(occ = c(1:4), #4 events
                        Phi = c(rep(moth.js.predicted$Phi$estimate, 3), NA),   # 3 survival estimates all the same
                        Nsuper = rep(moth.js.predicted$N$estimate + nrow(moth), 4), # Nsuper estimate + number of marked animals
                        pent = c(1-sum(moth.js.predicted$pent$estimate), moth.js.predicted$pent$estimate)) # Sum of all pent must be 1


N.derived <- data.frame(occ = c(1:4), #4 events
                        Phi = c(rep(moth.js.predicted$Phi$estimate, 3), NA),   # 3 survival estimates all the same
                        Nsuper = rep(moth.js.predicted$N$estimate + nrow(moth), 4), # Nsuper estimate + number of marked animals
                        pent = c(1-sum(moth.js.predicted$pent$estimate), moth.js.predicted$pent$estimate)) # Sum of all pent must be 1

N.derived$N <- NA
N.derived$N[1] <- (N.derived$Nsuper[1] * N.derived$pent[1])
for(i in 2:nrow(N.derived)){
  N.derived$N[i] <- (N.derived$N[i-1]*N.derived$Phi[i-1]) + (N.derived$Nsuper[i] * N.derived$pent[i])
}

N.derived

Step 3 (Partly. The code end on the line where I have problem)

detach("package:marked", unload=TRUE) # many of the function names are the same. unload `marked`
library(RMark)

moth.rmark.processed <- RMark::process.data(moth,
                                            model = "POPAN")

moth <- read.csv("moth.eps1.csv", colClasses = c(ch = "character"))
head(moth)

moth.rmark.processed <- RMark::process.data(moth,
                                            model = "POPAN")
Phi.time <- list(formula=~time)
p.dot <- list(formula=~1)

pent.sex <- list(formula=~sex)
N.sex <- list(formula=~sex)

library(RMark)
library(TMB)
moth.rmark <- mark(moth, model = "POPAN", model.parameters = list(Phi = Phi.time, p= p.dot, pent = pent.sex, N = N.sex),realvcv = TRUE)

PROBLEM: Last line code showing this error:

*Error in make.mark.model(data.proc, title = title, parameters = model.parameters, :

Error: Variable moth used in formula is not defined in data

Error in mark(moth, model = "POPAN", model.parameters = list(Phi = Phi.time, : Misspecification of model or internal error in code*

My best model is Phi(~time)p(~1)pent(~sex)N(~sex) . I assume that the problem is with sex variable because when I tried to change it to pent or N dot (~1) or time, the code was fine.

Thank you in advance for your kind help!