¡Descarga Scripts de R y más Apuntes en PDF de Biotecnología solo en Docsity!
COSAS DE LA ESTADISTICA DE PRIMERO
Introducció dades forma manual
Exemple de construcció d'una base de dades (data.frame) amb 4 variables, mitjançant
l'acumulació de 4 vectors en columna
Quan es resolguin els problemes concrets de la llista, caldrà:
1. afegir o suprimir columnes, segons necessitat
2. afegir o treure casos, segons necessitat
3. substituir els valors d'aquest exemple pels del problema
col1 <- c(1.1, 2.2, 3.3, 4.4, 5.5, 6, 7.7, 8, 9, 10) # dades exemple col 1, observeu la utilització del punt decimal i no la coma
col2 <- c(-1E3, -2E-2, -300E4, 4, -5, -6, -7, -8, NA, -10) # dades exemple col 2, amb exponents, la 9na dada missing
col3 <- rep(1,10) # dades exemple col 3, repetició 10 cops "1"
col4 <- 1:10 # dades exemple col 4, seqüència 1 a 10
en la instrucció següent afegir o treure els 'colx' que calgui
myData <- data.frame(col1, col2, col3, col4)
Opcional, si es desitja posar noms a les variables (les columnes)
colnames(myData) <- c("Variable_1", "Variable_2", "Repetició", "Seqüència")
Opcional, si es desitja posar noms als casos (les files)
rownames(myData) <- c("Cas 1", "Cas 2", paste("Cas",seq(3,10))) # observeu la generació dels noms 3 a 10
Introducció dades amb taula
Exemple de construcció d'un dataset important una taula Excel (convertida a csv)
Quan es resolguin els problemes concrets de la llista:
1. canviar el paràmetre 'myHeader' a FALSE si la 1ra fila de la taula csv no correspona
l'encapçalament
2. canviar el paràmetre 'myPath' per que apunti al directori on es troba el arxiu
3. canviar el paràmetre 'myFile' per especificar l'arxiu csv que conté les dades
4. canviar el paràmetre 'myNAspec' per que els missings/NA corresponguin a cel·les del csv
que continguin aquesta tira
paràmetres a modificar si s'escau
myHeader <- TRUE # posar FALSE si no hi ha fila d'encapçalament myPath <- '' # per defecte el directori on apunta RGui, si s'ha d'especificar seguir
l'exemple següent: 'C:/Documents and Settings/user/Mis
documentos/Exemples dead' myFile <- 'ExempleInicial.csv' # nom de l'arxiu, extensió inclosa myNAspec <- ' ' # així els NA corresponen a un espai en blanc
myData <- read.table(paste(myPath,myFile,sep=""),sep=";",header=myHeader, dec=",", na.strings=myNAspec);
Assumim que
1. el dataframe a treballar és myData (script 0.1 Entrada manual dades.r)
boxplot(myData[,myVars],main=myTitle,col=myColours)
histograma de la primera variable
hist(myData[,1],main='Histograma',col='palegreen')
scatterplot variables 1 i 4
plot(myData[,1],myData[,4],main='Diagrama de dispersió',col='red',xlab='una variable',ylab='una altra variable')
Proves amb 1 y 2 muestras (siguen una normal)
Assumim que
1. el dataframe a treballar és myData
2. test mitjana d'una mostra: paràmetres a especificar
myVar: índex de la variable a contrastar. Alternativament, es pot cridar també utilitzant
el nom de la variable
myH0: hipòtesi nul·la,
mySide: alternativa unil/bilat, confiança interval myConfLev
3. test mitjanes dues mostres:
myVars: índexs de les variables a estudiar (llista de 2 elements). Alternativament, noms
de les 2 variables
myDataPaired: dades independents/aparellades
myVarEq: assumim variables homogènies
mySide: alternativa unil/bilat, confiança interval myConfLev
4. test variàncies dues mostres:
myVars: indexs de les variables a estudiar (llista de 2 elements). Alternativament, noms
de les 2 variables
myH0: hipòtesi nul·la,
mySide: alternativa unil/bilat, confiança interval myConfLev
2. Prova per la mitjana d'una mostra normal
myVar <- 1 # 1ra variable del dataFrame myH0 <- 3 # testem si mu0 = 3 myConfLev <- 0.95 # interval del 95% mySide <- 'g' # H1 unilateral dreta (altres possibilitats 't'=bilateral, 'l'=unit. esquerra)
t.test(myData[,myVar],mu=myH0,alternative=mySide,conf.level=myConfLev)
o bé
t.test(myData$Variable_1,mu=myH0,alternative=mySide,conf.level=myConfLev)
PRUEBAS UNA DE Y DOS PROPORCIONES
1. test proporcions d'una mostra (asintòtic): paràmetres a especificar
myX: nombre de succesos
myN: nombre de repeticions
myH0: probabilitat sota la hipòtesi nul·la,
mySide: alternativa unil/bilat
myConfLev: confiança interval
2. test proporcions dues mostres (asintòtic):
myX: vector de dos posicions amb el nombre de succesos
myN: vector de dos posicions amb el nombre de repeticions
mySide: alternativa unil/bilat
myConfLev: confiança interval
3. test exacte una proporció:
myX: nombre de succesos
myN: nombre de repeticions
myH0: probabilitat sota la hipòtesi nul·la,
mySide: alternativa unil/bilat
myConfLev: confiança interval
Prova per una proporció (asintòtic)
myX <- 10 # 10 succesos myN <- 20 # 20 repeticions myH0 <- 0.4 # testem si p0 = 0. myConfLev <- 0.95 # interval del 95% mySide <- 'g' # H1 unilateral dreta (altres possibilitats 't'=bilateral, 'l'=unit. esquerra)
prop.test(x=myX,n=myN,p=myH0,alternative=mySide,conf.level=myConfLev)
Prova per a dues proporcions
myX <- c(10,15) # 10 i 15 succesos myN <- c(26,32) # 26 i 32 repeticions myConfLev <- 0.95 # interval del 95% mySide <- 't' # H1 bilateral (altres possibilitats 'g'=unit. dreta, 'l'=unit. esquerra)
Pruebas chi cuadrado
Assumim que
1. el dataframe a treballar és myData
2. test bondat ajustament: paràmetres a especificar
myVar: columna on son les freqüències
myProbs: vector de proporcions sota H
3. test homogeneïtat/independència: paràmetres a especificar
myVars: columnes corresponents a la variable categòrica de la taula de contingència
myProbs: vector de proporcions sota H
1. Prova bondat ajustament
.... dades exemple ....
myData <- data.frame(c(89,37,30,28,2))
.... paràmetres i crida mètode ....
myVar <- 1 # 1ra variable del dataFrame myProbs <- c(40,20,20,15,5) # vector probabilitats, es reescala tot a 1
chiTest <- chisq.test(myData[,myVar],correct = FALSE,p=myProbs,rescale.p= TRUE) print(chiTest$expected) print(chiTest)
2. Prova independència / homogeneïtat.
.... dades exemple ....
taula 4 files (variable categorica A) i 2 cols (variable categoria B)
col1 <- c(100, 200, 300, 400) # dades categoria B1: els valors son els comptatges de A1, A2 A3 i A4 en B col2 <- c(15, 25, 70, 10) # dades categoria B2: els valors son els comptatges de A1, A2 A3 i A4 en B myData <- data.frame(col1,col2) colnames(myData) <- c('B1','B2') rownames(myData) <- c('A1','A2','A3','A4')
.... paràmetres i crida mètode ....
myVars <- c(1,2) # 1ra i segona variables del dataFrame. Si B té k categories, posar c(1:k) chiTest <- chisq.test(myData[,myVars],correct = FALSE) print(chiTest$observed) print(chiTest$expected) print(chiTest)
y <- myData[,myVars[2]]
reg <- lm (y ~ x)
3. Diferents resultats
summary(reg) # resum general del ajust
anova(reg) # test coeficient de correlació reg$coefficients[1] # ordenada a l'origen reg$coefficients[2] # pendent
reg$residuals # residus
confint(reg,level=0.95) # interval de confiança dels coeficients. Canviar el 0.95 per altres confiances si s'escau
4. Prediccions
predict(reg,data.frame(x=c(7,-2)),interval='prediction') # cerquem la predicció dels valors x=7 i després x=- predict(reg,data.frame(x=c(7,-2)),interval='confidence') # cerquem la predicció dels valors x= i després x=-
5. Scatterplot i gràfics de diagnòstic
plot(x,y,main='Exemple scatterplot')
plot(reg)
6. Alternativa, si només interssa el test del coef. de correlació
cor.test(x,y)
AQUI EMPIEZA DEAD
ANOVA 1 FACTOR
Assumim que el dataframe a treballar és myData, amb dues columnes. La 1ra és el factor, la
segona la variable resposta
Modificar d'acord a les característiques del dataframe concret del problema (ubicacions de les
columnes i noms)
0.a Opcional: recordar que si s'escau es pot generar el vector amb els nivells del factor amb la
instrucció rep concatenada
per exemple, un factor amb quatre nivells que volem codificar amb 1, 2, 3 i 4 i que tenen
respectivament 4, 3, 5 i 6 rèpliques es genera amb
factor <- c(rep(1,4),rep(2,3),rep(3,5),rep(4,6))
millor que
factor <- c(1,1,1,1,2,2,2,3,3,3,3,3,4,4,4,4,4,4)
0.b Opcional: posar noms a les variables=columnes. Aqui agafem dos noms genèrics per fer
l'ANOVA
colnames(myData) <- c('factor', 'resp')
PREMISES DE LA ANOVA (CONDICIONS APLICACIO)
Comprovació opcional de les premisses de l'ANOVA
En aquest script assumim que:
1. el dataframe a treballar és myData, amb dues columnes.
2. la 1ra columna és el factor i porta com a nom 'factor', la segona la variable resposta, nom
columna 'resp'
3. s'ha calculat prèviament la taula anova i està emmagatzemada en l'objecte data.aov
Modificar d'acord a les característiques del dataframe concret del problema (noms de les
columnes), nom del anova etc
Opcional 1: 4 gràfics de control
plot(data.aov)
Opcional 2: test de homocedasticitat (igualtat de variàncies)
bartlett.test(resp ~ factor, myData)
Opcional 3: test de normalitat sobre els residus (NO sobre les dades originals!)
shapiro.test(data.aov$residuals)
ESTADISTICA DESCRIPTIVA ANOVA 1F
Estadistica descriptiva opcional Anova 1F
En aquest script assumim que:
1. el dataframe a treballar és myData, amb dues columnes.
2. la 1ra columna és el factor i porta com a nom 'factor', la segona la variable resposta, nom
columna 'resp'
Modificar d'acord a les característiques del dataframe concret del problema (nom i posició de
les columnes)
Opcional 1: Descriptiva variable resposta global i per grups.
summary(myData$factor) # rèpliques en cada nivell summary(myData$resp) # summary global by(myData$resp,myData$factor, summary) # summary per nivell by(myData$resp,myData$factor, sd) # desv. estandard per nivell
Opcional 2: boxplot segons nivell de factor
boxplot(resp ~ factor, myData)
Opcional 3: gràfic mitjanes vs variàncies
plot(tapply(myData$resp, myData$factor, var), tapply(myData$resp, myData$factor, mean), xlab='Variància resposta segons factor',
ANOVA 1 FACTOR ALEATORI
Assumim que el dataframe a treballar és myData, amb dues columnes. La 1ra és el factor, la
segona la variable resposta
Modificar d'acord a les característiques del dataframe concret del problema (ubicacions de les
columnes i noms)
0.a Opcional: recordar que si s'escau es pot generar el vector amb els nivells del factor amb la
instrucció rep concatenada
per exemple, un factor amb quatre nivells que volem codificar amb 1, 2, 3 i 4 i que tenen
respectivament 4, 3, 5 i 6 rèpliques es genera amb
factor <- c(rep(1,4),rep(2,3),rep(3,5),rep(4,6))
millor que
factor <- c(1,1,1,1,2,2,2,3,3,3,3,3,4,4,4,4,4,4)
0.b Opcional: posar noms a les variables=columnes. Aqui agafem dos noms genèrics per fer
l'ANOVA
colnames(myData) <- c('factor', 'resp')
0.c Obligatori: si no ha estat prèviament definit com a factor (en el sentit de R), convertir la
columna
myIndFactor <- 1 myData[,myIndFactor] <- as.factor(myData[,myIndFactor])
1. Càrrega del package nlme. Només cal fer-un cop per sessió
library(nlme)
2. sentències per fer l'ANOVA. Es guarda en l'objecte data.aov. Fins aqui, idèntic al anova 1F
d'efecte fix data.aov <- aov(resp ~ factor, myData) # generem l'objecte data.aov, un ANOVA 1F anova(data.aov)
3. si és significatiu, aleshores calcular les components de la variància
data.lme <- lme(resp ~ 1, random = ~ 1 | factor, data = myData,method = "REML") varcomp <- VarCorr(data.lme) varcompN <- as.double(varcomp)
4. opcional, calcular els percentatges de variabilitat
paste('Component variància factor', round(varcompN[1]/sum(varcompN[1:2])100,4),'%') paste('Component variància residual',round(varcompN[2]/sum(varcompN[1:2])100,4),'%')
TEST KRUSKAL WALLIS
