R: diferència entre les revisions

Salta a la navegació Salta a la cerca
4.863 bytes afegits ,  1 juny 2005
Cap resum de modificació
Línia 565: Línia 565:
A continuación se irán mostrando las distintas partes del código con una breve descripción de los procedimientos que realiza, para que cada cual vaya descubriendo por si sólo la funcionalidad y sintaxis de las distinas partes:
A continuación se irán mostrando las distintas partes del código con una breve descripción de los procedimientos que realiza, para que cada cual vaya descubriendo por si sólo la funcionalidad y sintaxis de las distinas partes:


1-
1-Definimos el nombre de la función, los parámetros que deberá contener y sus valores por defecto en caso que los haya:


interval<-function(model,digits=2,nivell=0.95,alfa=F) {
2-Creamos algunas variables internas a partir de los parámetros recibidos
nivell2<--round(qnorm((1-nivell)/2),digits=2)
co<-summary(model)$coef
tot<-c()
3-Comprobamos que estamos recibiendo un modelo de regresión válido:
if(class(model)[1]!="lm" & class(model)[1]!="glm" & class(model)[1]!="coxph")
{
ex.warn <- options()$warn
options(warn = 0)
warning("El objeto no es un modelo de regresion valido",call.=F)
options(warn = ex.warn)
}
4-Cálculos si hemos recibido un modelo de regresión lineal:
if(class(model)[1]=="lm")
{
df<-model$df
p<-(1-nivell)/2
punt<-round(co[,1],digits=digits)
inf<-round(confint(model,level=nivell)[,1],digits=digits)
sup<-round(confint(model,level=nivell)[,2],digits=digits)
tot<-cbind(punt,inf,sup)
dimnames(tot)<-list(row.names(summary(model)$coef),c("B","inf","sup"))
}
5-Cáculos si hemos recibido un modelo de riesgos proporcionales de Cox:
if(class(model)[1]=="coxph")
{
punt<-round(exp(co[,1]),digits=digits)
inf<-round(exp(co[,1]-nivell2*co[,3]),digits=digits)
sup<-round(exp(co[,1]+nivell2*co[,3]),digits=digits)
tot<-cbind(punt,inf,sup)
dimnames(tot)<-list(row.names(summary(model)$coef),c("HR","inf","sup"))
alfa<-T
}
6-Cálculos si hemos recibido un ''glm'' (logística, logbinomial, poisson):
if(class(model)[1]=="glm")
{
if(summary(model)$family[1]=="gaussian")
{
ex.warn <- options()$warn
options(warn = 0)
warning("Utiliza la funcion lm en lugar de glm para los modelos de regresion lineal",call.=F)
options(warn = ex.warn)
}
else
{
punt<-round(exp(co[,1]),digits=digits)
inf<-round(exp(co[,1]-nivell2*co[,2]),digits=digits)
sup<-round(exp(co[,1]+nivell2*co[,2]),digits=digits)
tot<-cbind(punt,inf,sup)
if(summary(model)$family[1]=="binomial" & summary(model)$family[2]=="logit")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("OR","inf","sup"))
}
if(summary(model)$family[1]=="binomial" & summary(model)$family[2]=="log")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("RR","inf","sup"))
}
if(summary(model)$family[1]=="poisson")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("IRR","inf","sup"))
}
}
}
7-Preparamos e imprimimos los resultados:
if(!is.null(tot) & alfa==T)
{
return(tot)
}
if(!is.null(tot) & alfa==F)
{
tot<-tot[2:dim(tot)[1],]
if(length(tot)>3)
{
return(tot)
}
else
{
tot<-t(tot)
row.names(tot)<-row.names(summary(model)$coef)[2]
return(tot)
}
}
}
A continuación se presenta el código completo. Su edición con [[Emacs]] o [http://www.winedt.com WinEdt] facilita su comprensión al colorear las funciones y los parentésis y claudators.
interval<-function(model,digits=2,nivell=0.95,alfa=F) {
nivell2<--round(qnorm((1-nivell)/2),digits=2)
co<-summary(model)$coef
tot<-c()
if(class(model)[1]!="lm" & class(model)[1]!="glm" & class(model)[1]!="coxph")
{
ex.warn <- options()$warn
options(warn = 0)
warning("El objeto no es un modelo de regresion valido",call.=F)
options(warn = ex.warn)
}
if(class(model)[1]=="lm")
{
df<-model$df
p<-(1-nivell)/2
punt<-round(co[,1],digits=digits)
inf<-round(confint(model,level=nivell)[,1],digits=digits)
sup<-round(confint(model,level=nivell)[,2],digits=digits)
tot<-cbind(punt,inf,sup)
dimnames(tot)<-list(row.names(summary(model)$coef),c("B","inf","sup"))
}
if(class(model)[1]=="coxph")
{
punt<-round(exp(co[,1]),digits=digits)
inf<-round(exp(co[,1]-nivell2*co[,3]),digits=digits)
sup<-round(exp(co[,1]+nivell2*co[,3]),digits=digits)
tot<-cbind(punt,inf,sup)
dimnames(tot)<-list(row.names(summary(model)$coef),c("HR","inf","sup"))
alfa<-T
}
if(class(model)[1]=="glm")
{
if(summary(model)$family[1]=="gaussian")
{
ex.warn <- options()$warn
options(warn = 0)
warning("Utiliza la funcion lm en lugar de glm para los modelos de regresion lineal",call.=F)
options(warn = ex.warn)
}
else
{
punt<-round(exp(co[,1]),digits=digits)
inf<-round(exp(co[,1]-nivell2*co[,2]),digits=digits)
sup<-round(exp(co[,1]+nivell2*co[,2]),digits=digits)
tot<-cbind(punt,inf,sup)
if(summary(model)$family[1]=="binomial" & summary(model)$family[2]=="logit")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("OR","inf","sup"))
}
if(summary(model)$family[1]=="binomial" & summary(model)$family[2]=="log")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("RR","inf","sup"))
}
if(summary(model)$family[1]=="poisson")
{
dimnames(tot)<-list(row.names(summary(model)$coef),c("IRR","inf","sup"))
}
}
}
if(!is.null(tot) & alfa==T)
{
return(tot)
}
if(!is.null(tot) & alfa==F)
{
tot<-tot[2:dim(tot)[1],]
if(length(tot)>3)
{
return(tot)
}
else
{
tot<-t(tot)
row.names(tot)<-row.names(summary(model)$coef)[2]
return(tot)
}
}
}


=== Documentación de funciones ===
=== Documentación de funciones ===

Menú de navegació