Hi, I am relatively new to RMark but have been running my data using MARK for a bit longer. I'm modeling 4 primary occasions with various secondary lengths. We are trapping prairie dogs essentially for two primary periods a year. We have finished two field seasons (4 primary sessions) and I have one more to go). I want to move it into RMark to make it easier to add new data to models and run heterogeneity models as my design matrices were getting very complicated to produce. I am interested in looking at survival rates but I also need abundance estimates to report to the parks I am trapping on. I've been attempting to write my code by fudging around a friends code from a previous project but I seem to be running into a couple of problems:
1) For now I want to just model no movement (when I feel more comfortable I'll probably add in random movement which I had back in MARK but I'm trying to keep things as simple as I can for now), but with the code I am currently running, I am still getting Gamma' to come up .5 while Gamma'' is 0. What am I doing wrong?
2) Am I modeling p=c correctly? I think I read from others coding on here this is how I should be doing it but I just wanted to double check? We only catch about 80-120 animals each primary occasion so our p's tend to be low and as a result I need to keep p=c if I'm going to model time for sake of having too many parameters.
I am attaching the code that I am working with and any help would be greatly appreciated. Thanks so much,
Amanda
#RMark analysis
#Load RMark package
library(RMark);
rm(list=ls())
cleanup(ask=FALSE)
#inport .inp
bentsadults=convert.inp("C:/R_robust/inpfiles/bentsadultsrobust.inp",
group.df=data.frame(sex=c("Male", "Female")))
bentsadults.robust.run=function()
{
#Process the data so models can be analyzed
#Time intervals for 4 primary occasions (300+ collectors)
#use the format from MARK where you put in the x/60 for time intervals
time=c(0,0,0,0,0,0,.733,0,0,0,0,4.23,0,0,0,0,0,0,0,0,0,0,0,0,1.47,0,0,0,0,0,0,0,0,0,0,0,0,0)
bentsadults.proc=process.data(bentsadults,model="RDFullHet",time.intervals=time,
groups="sex")
summary(bentsadults.proc)
#Create the design data
bentsadults.ddl=make.design.data(bentsadults.proc)
#Create parameter specifications for formulas
#S.dot=list(formula=~1)
#S.group=list(formula=~group)
S.time=list(formula=~time)
#S.grpbytime=list(formula=~group*time)
#S.grpplustime=list(formula=~group+time)
#GammaDoublePrime.0=list(formula=~1,fixed=0)
#GammaPrime.0=list(formula=~1,fixed=0)
GammaDoublePrime.dot=list(formula=~1,fixed=0,share=TRUE)
#p.dot=list(formula=~1)
#p.time=list(formula=~time)
#p.group=list(formula=~group)
#p.time.group=list(formula=~time*group)
#p.timeplusgroup=list(formula=~time+group)
#p.dotsession=list(formula=~session)
p.dotc=list(formula=~session,share=TRUE)
p.dotgroupc=list(formula=~session:group,share=TRUE)
p.dottimec=list(formula=~session:time,share=TRUE)
p.timegroupc=list(formula=~session:time*group,share=TRUE)
p.timeplusgroupc=list(formula=~session:time+group,share=TRUE)
#c.dot=list(formula=~1)
#c.time=list(formula=~time)
#c.group=list(formula=~group)
#c.time.group=list(formula=~time*group)
#c.timeplusgroup=list(formula=~time+group)
#c.dotsession=list(formula=~session)
pi.dot=list(formula=~1)
pi.group=list(formula=~group)
pi.session=list(formula=~session)
pi.none=list(formula=~1, fixed=0)
#N.dot=list(formula=~1)
#N.group=list(formula=~group)
N.session=list(formula=~session)
N.grpbysession=list(formula=~group*session)
# create model list
cml=create.model.list("RDFullHet")
#row.names(cml)=c(1:length(cml$p))
# run and return models
results=mark.wrapper(cml,data=bentsadults.proc,ddl=bentsadults.ddl,adjust=TRUE)
return(results)
}
robust.results=bentsadults.robust.run()
robust.results
results.table=robust.results$model.table