El modelo de CJS (Cormack-Jolly-Seber) es el básis de muchos de los modelos de CMR para poblaciones abiertas que vamos a tratar, por lo que un firme entendimiento de este modelo y sus aplicaciones es esencial. Voy a ilustrar el modelo de CJS con un ejemplo bien conocido de datos de anillar y recapturar en vivo para el Dipper Europea. Este ejemplo de datos ha sido utilizado muchas veces para ilustrar CJS y era un ejemplo clave en Lebreton y col. (1992) Monografía Ecológica. El conjunto de datos consistió de 7 ocasiones de marcado y recaptura para Dippers machos y hembras y se utilizará para ilustrar:
La estructura básica del modelos (vs. cohorte) tiempo específico de las probabilidades de supervivencia y de recaptura
El uso de covariables de agrupación (en este caso, el sexo), un caso especial de covariables individuales
El uso de covariables tiempo específicos para modelar la variación en la supervivencia o recuperar
Evaluación de ajuste de los modelos y comparación de modelos alternativos
Les llevaré a través de los pasos básicos de la construcción de modelos de CJS en RMark continuación. Más detalles se proporcionan en el Apéndice C de línea introducción agradable de Evan Cooch a MARK.
Los datos de los Dippers
En primer lugar, asegúrese de que la biblioteca de RMark se ha cargado correctamente en la memoria. Si no está ya cargado, utilice el comando
> require (RMark)
Los datos del cazo se almacenan en una trama de datos de RMark. Recordemos que hemos comentado anteriormente cómo dar formato a los datos para MARK o RMark, y voy a asumir que sus datos están en este formato, al igual que los datos Dipper. Cargue estos datos utilizando el comando (ver el script de R para conveniencia)
> data (dipper)
Podemos mostrar los datos
> dipper
ch sex
1 0000001 Female
2 0000001 Female
3 0000001 Female
4 0000001 Female
5 0000001 Female
6 0000001 Female
7 0000001 Female
8 0000001 Female
......
263 0110000 Male
264 0110110 Female
265 0111000 Female
266 0111000 Male
267 0111100 Female
268 0111100 Female
etc
Lo que vemos es que cada línea contiene una cadena de caracteres de la historia de captura (de 1 y 0) (7 columnas en total), seguido por un carácter "Male" o "Female" que denota el sexo del ave individual. Las ocasiones de captura están atrapando anual esfuerzos a partir de 1980. Utilizamos esta información junto con el conocimiento de que vamos a utilizar el modelado CJS para "procesar" los datos en un formulario que proporcionará información más significativa para el análisis.:
#process dipper data to CJS format
dipper.processed=process.data(dipper,model="CJS",begin.time=1980,groups="sex")
Esto, por ejemplo, nos permitirá usar el "sexo" en lugar de una variable de grupo genérico para el modelado, y hará un seguimiento de las ocasiones por años naturales reales, por lo que la salida más interpretable. Por último, vamos a crear un "datos de diseño" (*. ddl) por el comando
>
> dipper.ddl=make.design.data(dipper.processed)
Aunque no es absolutamente necesario (un archivo ddl se crea de forma predeterminada en la marca () de comandos), este paso es una buena práctica y nos dará algunos detalles útiles que podemos dejar para más tarde.
Construir modelos
Los modelos son (aparentemente) fácil de crear en RMark; podemos crear nuestro primer modelo sencillo con sólo usar las anteriores 2 ingredientes (los datos procesados y DDL) como argumentos en la marca () del sistema:
> my_first_model<-mark(dipper.processed,dipper.ddl)
Varias cosas suceden como resultado de este comando:
Usted recibirá una buena cantidad de salida de pantalla que describe el modelo, las estimaciones, y así sucesivamente.
Los resultados (y datos) todos serán almacenados en un objeto llamado my_first_model
Varios archivos con prefijos mark001, mark002, etc se crearán en el directorio de trabajo
Vamos a hablar sobre todo acerca de los 2 primeros de éstos; el tercero será útil para las aplicaciones más avanzadas (o si se quiere portar su análisis a la interfaz gráfica de usuario de MARK).
Veamos primero a los primeras líneas de los resultados:
> my_first_model<-mark(dipper.processed,dipper.ddl)
STOP NORMAL EXIT
Output summary for CJS model
Name : Phi(~1)p(~1)
Npar : 2
-2lnL: 666.8377
AICc : 670.866
Estas líneas indican, en primer lugar, que se logró una ejecución normal (sin errores); segundo, nos dicen qué modelo específico fue estimado en análisis CJS, y luego nos proporcionan una serie de estadísticas adicionales (número de parámetros, la desviación, el AIC). En particular, note que el modelo
Phi (~ 1) p (~ 1)
indica que sólo las intercepciones se han utilizado para modelar la variación en Phi y p. Usted puede preguntarse "¿de dónde viene este modelo" y la respuesta es que "es el modelo por defecto utilizado por RMark¨ si no se especifica una estructura modelo CJS (que hizo en el paso de proceso). Como regla general, recomiendo la creación explícita de los modelos, que vamos a hacer a continuación. Pero primero vamos a examinar a algunos de la otra salida. Primero, tenemos los estimados "beta". Estas son los estimados de los parámetros en la escala logit, y en este caso son simplemente los 2 intercepciones (por Phi y p), junto con las SE y los intervalos de confianza.
Beta
estimate se lcl ucl
Phi:(Intercept) 0.2421484 0.1020127 0.0422034 0.4420933
p:(Intercept) 2.2262661 0.3251094 1.5890517 2.8634805
Después de esta son una serie de tablas que contienen los estimados de los parámetros "reales" (es decir, hacia atrás-transformado) de Phi y p. Observe que estos son en la forma de matrices diagonales, con las filas y columnas marcadas por año (1980-1986). Específicamente, las filas representan la cohorte de la suelata (año en el que un grupo de animales marcados fueron soltados juntos y las columnas, ya sea para las ocasiones de encuentro (p) o el intervalo en el que se produjo la supervivencia (Phi). Hay una matriz separada para cada tipo de parámetro, y para cada grupo (el sexo), por lo que 4 en total (2 de Phi y 2 para p) en este caso. los valores de estas matrices se alinean con el índice del parámetro correspondiente en lo que se conoce como las matrices de parámetros (índice PIMS). vamos a entrar en esto en más detalle a continuación, pero brevemente que representaremos cómo los parámetros son potencialmente variando los grupos, ocasiones, o cohortes (o combinaciones de estas cosas).
Una cosa que ustedes pueden haber notado es que todos los valores de Phi y p en las matrices son idéntico. Esto es debido a que el modelo -. que vamos a especificar en la Matriz de Diseño -. especifica que no hay ninguna variación en alguno de estos factores. RMark sin embargo parte de una premisa muy general de que es posible que desee permitir a todos tipo de variación en estos parámetros - por lo que estos son marcadores de posición para permitir que la variación. De hecho, una versión simplificada de la PIMS tendría sólo un único valor para cada uno de Phi y p (2 parámetros), y sería proporcionar los mismos resultados.
Vamos a crear algunos modelos adicionales, más interesantes. Estos modelos se basan en la forma, dada la estructura de datos, podríamos esperar que los parámetros que varían con ocasiones (tiempo) o entre grupos (sexo). Por ahora, no vamos a permitir que la variación debida a la cohorte de la suelta, es decir, vamos a utilizar supuestos CJS estándar que se capturan animales sobrevivientes de cohortes anteriores y sobrevivir con las mismas probabilidades como las nuevas versiones. Más tarde nos relajaremos este bajo supuestos que permiten la cohorte (por lo general esto significa años de edad) de la dependencia.
Vamos a representar este tipo de variación de cada parámetro mediante la creación de fórmulas que será familiar por ahora - porque siguen el formato R estándar para modelos lineales generalizados. Para Phi tenemos:
> Phi.dot=list(formula=~1)
> Phi.sex=list(formula=~sex)
> Phi.t=list(formula=~time)
> Phi.sex.t=list(formula=~sex+time)
> Phi.sex.T=list(formula=~sex*time)
representando respectivamente sin variación (sólo en el origen), la variación por sexo solamente, sexo + modelo aditivo tiempo, el sexo * modelo interactivo del tiempo. Asimismo para p tenemos
> p.dot=list(formula=~1)
> p.sex=list(formula=~sex)
> p.t=list(formula=~time)
> p.sex.t=list(formula=~sex+time)
> p.sex.T=list(formula=~sex*time)
Podemos poner esto juntos para formar 25 modelos alternativos (5 Phi * 5 p). Sólo voy a mostrar los comandos de los primeros 5, pero el resto están en el archivo RSCRIPT creé.
>mod1=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.dot))
>mod2=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.sex))
>mod3=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.t))
>mod4=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.sex.t))
>mod5=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.sex.T))
Después de ejecutar todos los 25 modelos, podemos usar la función incorporada collect.models() para solucionarlos por los valores de la AIC. Esta función sólo se ve en todos los modelos que se han ejecutado y que están todavía en la memoria, así que recuerde que en realidad corrió mod1 dos veces (la primera vez como "mymodel" con los valores de intersección por defecto). No queremos que este aparezca dos veces en las mesas de la AIC, así que vamos a quitar uno (no importa cuál)
>rm(myexample)
>collect.models()
He salvado la salida de la pantalla a un archivo de texto. No haga caso de la primera columna (un índice arbitrario para cada modelo) y el enfoque de los valores de la AIC y pesos de la AIC. Los modelos están clasificadas por valores ascendentes y descendentes de los pesos del AIC, con los modelos de "mejores" en primer lugar. Es evidente que el modelo simple Phi (~ 1) p (~ 1) clasifica mejor, pero es seguido de cerca por los modelos que permiten la variación sexual en Phi o p. Por otro lado, la mayoría de los modelos que incorporan efectos de tiempo se orden bastante bajo. Discutiremos más adelante cómo hacer frente a la incertidumbre del modelo (por ejemplo, la selección del modelo, el modelo de estimación promedio).
Tenga en cuenta que todavía tenemos el modelo de objetos para cualquier modelo que no borramos (quitar), y estos se pueden utilizar para generar resúmenes y otra información. por ejemplo
>summary(mod10)
produce una tabla de valores de resumen para mod10
> mod10 $ results
produce un conjunto detallado de los resultados, incluso las matrices de varianza-covarianza.
> mod10 $ output
nos dice donde están los archivos MARK están asociados a este modelo (muy útil si queremos convertir el anális a MARK
>attributes(mod10)
da un conjunto completo de atributos de este objeto.
Ahora, algunos de ustedes han utilizado MARK antes y pueden incluso haber ejecutar el ejemplo Dipper. Si es así, usted sabe que usted podría utilizar la GUI de MARK para crear los modelos anteriores (y más) como "modelos predefinidos", para manipular las cartas de PIM, las matrices de diseño, etc. Así que tal vez usted no está impresionado de que nos acabamos de escribir un montón de script para hacer lo mismo. En mi mente, hay varias razones por las que podría en realidad prefieren el enfoque RMark:
Más o menos RMark se obliga enfrentarse con lo que el modelo es, en notación de que está bastante claro, algo que puede no suceder en absoluto si usted auto-mágicamente crear un montón de modelos
A medida que crea cada vez modelos más complejos, y muchos de ellos, es muy probable que llegar a un punto en el que la creación de modelos en la interfaz de MARK vuelve tedioso y propenso a errores.
La interfaz de usuario gráfica (GUI) en MARK pone una cierta cantidad de sobrecarga adicional sobre los recursos de computación y en casos extremos puede conducir a falla del sistema. RMark elimina este overhead al interactuar directamente con los mecánicos de MARK.
Una vez que se han construido un conjunto básico de las declaraciones de modelos modelo, muchos otros modelos se pueden crear de forma rápida mediante la combinación de secuencias de comandos, como se ilustra en el ejemplo anterior. Una vez construí las 5 afirmaciones que describen la variación en Phi y p, y se incorporan un conjunto de éstos en función de mark(), que era capaz de crear rápidamente 25 modelos por copiar, pegar, y una cantidad relativamente modesta de edición de texto. Me tomó solo 10 minutos, y yo escribo lentamente.
La escritura de secuencias de comandos construye el carácter y te hará una mejor persona.
PIMS y matrices de diseño: una mirada más cercana
Como se mencionó anteriormente, RMark permite el modelado extremadamente flexible y general de CJS y otros datos. En consecuencia, hay un cierto grado de complejidad incrustada que el usuario medio puede o no puede necesitar saber. Sin embargo, hay algunas cosas básicas que usted debe conocer. Voy a tratar de presentarlos como unconfusing una manera de lo posible, volviendo al ejemplo del cazo para ilustrar varios puntos. Tomando este ejemplo, vamos a considerar lo que se logra mediante los siguientes pasos
> data_dir<-"C:/users/mike/Dropbox/teaching/WILD8390/spring2014/CJS"
> setwd(data_dir)
> rm(list=ls())
>
> library(RMark)
> #load dataframe containing dipper example
> data(dipper)
> #process dipper data to CJS format
Las primeras líneas se hicieron cargo de algunas tareas básicas, tales como la definición de un directorio de trabajo, la carga de la biblioteca RMark, y cargar los datos de Dipper. La línea
> dipper.processed=process.data(dipper,model="CJS",begin.time=1980,groups="sex")
establece que estamos siguiendo una estructura del modelo de CJS básico (frente a un otro), que la primera ocasión el tiempo es ser etiquetado como 1980, y que hay una variable de agrupación llamado "sexo". La línea
> dipper.ddl=make.design.data(dipper.processed)
crea una estructura de datos de diseño muy general para estos datos y este modelo. Si examina esta estructura de datos (capturado aquí en un archivo de texto) se puede ver que 1) se compone de 2 sublistas, una para Phi y otro para p, y 2) para cada uno hay 42 filas. Las filas corresponden a la "matriz no simplificada de indices de parámetros (PIM)" que veremos en detalle más adelante; las columnas son diferentes de los parámetros y el modelo de índices y factores que potencialmente permiten el modelado por grupo, tiempo, cohorte, la edad y otros factores. Esta estructura de datos se crea utilizando la función make.design.data () , y se puede agregar a la tarde, como veremos (por ejemplo, para permitir la inclusión de otras variables). Cuando construimos modelos que incluyen varios términos (tiempo, sexo, etc) estos términos deben hacer referencia a las columnas de esta trama de datos (por ejemplo, si usted mal hechizo sex como Sex Rmark no sabrá lo que está hablando).
Vamos a construir unos modelos para mostrar cómo funciona esto en la práctica. Una vez más, podemos definir un puñado de fórmulas para Phi y p
> Phi.sex=list(formula=~sex)
> Phi.t=list(formula=~time)
> Phi.sex.t=list(formula=~sex+time)
> p.sex=list(formula=~sex)
> p.t=list(formula=~time)
> p.sex.t=list(formula=~sex+time)
y utilizarlos para construir algunos modelos. Vamos a construir 3 modelos, en orden de complejidad decreciente:
>#create a model with sex+time for both phi and p
>mod1=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.sex.t,p=p.sex.t))
>#create a simpler model sex only for phi sex+time for p
>mod2=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.sex,p=p.sex.t))
#create an even model with phi.sex and p.sex
>mod3=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.sex,p=p.sex))
> # crear un modelo con el sexo + tiempo para ambos phi yp
> mod1 = marca (dipper.processed, dipper.ddl, model.parameters = list (Phi = Phi.sex.t, p = p.sex.t))
> # crear un modelo más simple de sexo sólo para phi sexo + tiempo para p
> mod2 = marca (dipper.processed, dipper.ddl, model.parameters = list (Phi = Phi.sex, p = p.sex.t))
# crea un modelo aún con phi.sex y p.sex
> mod3 = marca (dipper.processed, dipper.ddl, model.parameters = list (Phi = Phi.sex, p = p.sex))
Estos contienen 84 índices de parámetros, lo que permite la especificidad de grupo, el tiempo, y de cohortes, tanto en Phi y p, mucho más detalle de lo que va a ser apoyado por los datos o incluso estadísticamente identificable. De hecho, las PIMs no simplificadas para todos nuestros modelos tendrán esta misma dimensión. Pero, incluso nuestros modelos más complejos se van a simplificar esta estructura algo por el colapso de algunas dimensiones. Una simplificación de serie en los modelos CJS es especificar que Phi y p dependen del tiempo (ocasión), pero no cohorte de la suelta. Podemos mostrar los PIMS simplificados por
>PIMS(mod1,"Phi",simplified=T)
>PIMS(mod1,"p",simplified=T)
Al hacer esto, verá que los índices cambian a través de las columnas (ocasiones) en cada grupo, pero ya no a través de las filas; esto se traduce en una reducción de 84 a 24 índices de los parámetros, y representa el grupo completo * modelo CJS tiempo.
Sin embargo, no hemos terminado - el diseño opera más en el PIMS para crear el modelo especificado en el script del modelo. Recordemos que este modelo es el sexo + tiempo para ambos Phi y p, y por lo tanto carece de la interacción plazo. Podemos ver el diseño
>#rows for Phi
>mod1$design.matrix[1:12,]
>#rows for p
>mod1$design.matrix[13:24,]
Usted verá que esta matriz de diseño cuenta con 24 filas (1 para cada parámetro en los PIM), pero sólo 14 columnas. Las columnas se refieren a los parámetros del modelo lineal (logit-escala) utilizado para crear el modelo. Por ejemplo, todas las filas de la phi de tener un 1 en la columna de la intersección, un 1 o un 0 en la columna de sexo (0 si es mujer, 1 si es hombre), y un 1 o un 0 en cada vez ocasión que denota la columna cuya ocasión el tiempo se está modelando (la primera ocasión en la que represente la intercepción).
Vamos a continuar con el modelo un poco más simple, mod2. Este modelo considera el sexo, pero la variación no es tiempo para Phi, y el sexo + tiempo de variación de p. Visualizar la estructura PIM simplificado para este modelo
>#simplified PIM structure
>PIMS(mod2,"Phi",simplified=T)
>PIMS(mod2,"p",simplified=T)
Usted verá que sólo los 2 índices de parámetros se necesitan para Phi (uno para cada sexo), pero 12 son necesarios para p, para 14 en total. La matriz de diseño entonces funciona en este modelo más simple; tiene 14 filas (una para cada índice de parámetro) pero sólo 9 columnas: una intercepción y columna sexo por Phi; interceptar, sexo y tiempo de columnas para p. Si usted recibe un resumen de las estimaciones de este modelo
>summary(mod2)
verás que hay 9 parámetros "Beta" estimados; estas estimaciones corresponden a la intersección logit escala y coeficientes para cada una de las columnas de la matriz de diseño, y se utilizan para derivados de las estimaciones de los parámetros "reales" por tranformation logit inversa. Por último, tenga en cuenta mod3, con efectos de grupo sólo para Phi y p. La estructura de PIM simplificado
>PIMS(mod3,"Phi",simplified=T)
>PIMS(mod3,"p",simplified=T)
Sólo cuenta con 4 índices, y la matriz de diseño
> mod3$design.matrix
tiene 4 filas y 4 columnas, por lo que, de hecho, la matriz de diseño no ha cambiado nada (hay 4 parámetros estimados).
Todo esto parece bastante complicado, pero en realidad es relativamente sencillo y se realiza detrás de las escenas para usted si usted sigue estos pasos:
Especificar correctamente la estructura de datos y el modelo de la cuenta de proceso
Asegúrese de que todos los términos del modelo que se utilizan se incluyen en los datos de diseño (de nuevo, esto se hace de forma automática para los casos simples, pero requerirá una cierta manipulación de los más complejos)
Escriba los términos correctamente para los modelos logit lineal
Hemos de tener en cuenta que
El modelo real es el que resulta de los cálculos en las matrices de diseño-PIM
El número de parámetros estimables en el modelo es siempre el mismo que el número de columnas de la matriz de diseño.
He incluido los comandos anteriores en un archivo de script de R para su conveniencia.
Incluyendo covariables en el modelo
Hay 2 tipos diferentes de covariables que podemos considerar en un análisis de CJS y se manejan de manera diferente.
Covariables diseño se relacionan con la estructura del modelo especificado y incluyen el tiempo, cohortes y grupos (por ejemplo, el sexo en el ejemplo del cazo). Estos son manejados por la manipulación de la estructura de datos de diseño (ddl)
Covariables individuales son atributos que varían de un animal individual a los animales (por ejemplo, la masa después de su captura). Covariables individuales, permanecen en el archivo histórico encuentro individuo (como columnas adicionales en ese archivo)
Como se notó en el Apéndice C, esto no es una separación verdaderamente clara, ya que, obviamente, algunos atributos como el sexo -varían de un animal a otro. Para la mayoría de los propósitos, las covariables de diseño son atributos que son categóricos y que se prestan a la agrupación de los animales, y las covariables individuales son atributos numéricos (longitud, masa, etc.). Aquí, he utilizado el ejemplo de los Dippers para ilustrar un par de maneras que covariables diseño pueden ser incluidos.
Lebreton y col. (1992) incorporó uno de los primeros (si no el primero) usos de una covariable tiempo específico en el análisis CJS. Observaron que los años entre 1980 y 1986 eran o años de las inundaciones (1981, 1982) o los años normales, y crearon un indicador variable = 1 si el año era de inundaciones y 0 si era normal. El código de incorporar esto en la dll requiere algunos datos bastante simples R trucos de manipulación:
> dipper.ddl$Phi$Flood=0
> dipper.ddl$Phi$Flood[dipper.ddl$Phi$time==1981 | dipper.ddl$Phi$time==1982]=1
La dll modificado se utiliza luego de una serie de fórmulas modelo y definiciones para crear y comparar 3 modelos alternativos;
>Phi.Flood=list(formula=~Flood)
>Phi.dot=list(formula=~1)
>Phi.time=list(formula=~time)
>p.dot=list(formula=~1)
>dipper.phi.flood.p.dot=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.Flood,p=p.dot))
>dipper.phi.t.p.dot=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.time,p=p.dot))
>dipper.phi.dot.p.dot=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.dot))
>collect.models()
Los resultados indicaron pruebas bastante fuerte (un peso de 0,89 AIC) para inundación como un factor de predicción de la supervivencia, con este modelo de realizar mejor que sea el modelo Phi.dot o el modelo genérico con variación en el tiempo (Phi.time). El estimdao de beta -0.5599740 (IC 95% = -0.9840706, -0.1358773) indica que la supervivencia fue menor en las inundaciones que en años normales.
También usé la variable artificial "esfuerzo", construido en el Apéndice C para ilustrar cómo algo se pueden incorporar covariables más complicados. Aquí una variable esfuerzo de captura fue creado en una trama de datos
> df=data.frame(group=c(rep("Female",6),rep("Male",6)),time=rep(c(1981:1986),2), effort=rep(c(10,5,2,8,1,2),2))
>df
group time effort
1 Female 1981 10
2 Female 1982 5
3 Female 1983 2
4 Female 1984 8
5 Female 1985 1
6 Female 1986 2
7 Male 1981 10
8 Male 1982 5
9 Malco 1983 2
10 Male 1984 8
11 Male 1985 1
12 Male 1986 2
Esto luego se fusionó con la dll utilizando un built-in merge_design.covariates función ():
> Dipper.ddl $ p = merge_design.covariates (dipper.ddl $ p, df, bygroup = TRUE)
Entonces fórmulas de parámetros se construyen que implica este covariable como un predictor de p, y se utilizan para construir algunos modelos adicionales
> p.effort.plus.sex=list(formula=~effort+sex)
> p.effort=list(formula=~effort)
>
> dipper.phi.dot.p.effort.sex=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.effort.plus.sex))
>dipper.phi.dot.p.effort=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.dot,p=p.effort))
>dipper.phi.flood.p.effort=mark(dipper.processed,dipper.ddl,model.parameters=list(Phi=Phi.Flood,p=p.effort))
>collect.models()
En este ejemplo, el modelo Phi(~ Flood) p (~ esfuerzo) ahora termina siendo uno de los modelos de gama alta, lo que indica algún tipo de apoyo para la covariable esfuerzo.
He incluido estos ejemplos en un script de R por conveniencia.
Bondad de ajuste,selección de modelos, y inferencia de modelos múltiples
Bondad de ajuste
Como se ha visto en los ejemplos anteriores hay una serie de aproximaciones para la bondad de ajuste. Modelos de CJS en realidad han recibido mucha atención en esta área, y algunos procedimientos bien desarrollados están disponibles. Ya he mencionado la aproximació de la desviación / df, y que puedo trabajar razonablemente bien para los modelos CJS.
El programa RELEASE lleva a cabo una serie de pruebas de tabla de contingencia de chi-cuadrado para los modelos CJS estándar, en los cuales se produce una bondad de ajuste de prueba para los grupos generales * modelo de tiempo (el modelo más complejo que hemos creado hasta ahora). Los comandos para llevar a cabo la prueba de gof RELEASE en RMark para el ejemplo son:
> dipper.processed=process.data(dipper,groups=("sex"))
> release.gof(dipper.processed)
RELEASE NORMAL TERMINATION
Chi.square df P
TEST2 7.5342 6 0.2743
TEST3 10.7735 15 0.7685
Total 18.3077 21 0.6295
La prueba sugiere muy buen ajuste del modelo. En comparación un cálculo de la desviación / df para ese modelo (mod25) proporciona
> mod25$results$deviance/mod25$results$deviance.df
[1] 3.761788
lo que sugiere cierta falta de ajuste o sobredispersión. Vamos a continuar esta labor con una prueba de ajuste de bootstrap paramétrico. Para implementar se necesita primero cargar un objecto de R que contiente 3 functions necesario para el bootstrap. Guardar este objeto en su directorio. El script de R contiene los comandos necesario para correr el boostrap y las pruebas de RELEASE y desviación / df. Los resultados sugerien la necesidad de un ajuste de c-hat mas pequeño que desviación / df .
Inferencia de modelos múltiples
En general, no vamos a querer asumir que un modelo único es operativo, sino más bien hacer inferencias (parámetros de estimación, hacer comparaciones, hacer predicciones), teniendo en cuenta la incertidumbre del modelo según lo expresado por los pesos de la AIC. Ilustro con un ejemplo de Parus major, un carbonero anillado (congregado) en España como parte de un estudio de la selección para "Tamaño de lazo" (tamaño de una banda de pecho negro), una marca en los hombres que parece que se correlaciona con la supervivencia (Senar et al., en prensa). Aquí están los datos (en formato MARK) para adultos y subadultos carboneros, y el código R para dar formato a los datos para RMark, lea covariables individuales, conducta análisis CJS por edad, y calcular y trazar las predicciones del modelo-un promedio de supervivencia.
Siguiente- Modelos de edades múltiples