Creación de funciones
cubos <- function(x) x**3
#Y la usamos o llamamos como
cubos (1:5)
## [1] 1 8 27 64 125
## vector
v1 <- seq(from = 0, to= 20,by= 2); v1
## [1] 0 2 4 6 8 10 12 14 16 18 20
# Funcion en el vector
cubos(v1)
## [1] 0 8 64 216 512 1000 1728 2744 4096 5832 8000
Estandarización de datos
estandarizar <- function(x) {
med= mean(x)
dv= x-med # desvios
z= dv/sd(x)
par(mfrow=c(2,2))
hist(x)
hist(dv)
hist(z)
return(list(med.data=med,m.dv=mean(dv),var.x=var(x),var.dv=var(dv),
m.z=mean(z), var.z=var(z),z=z))
}
set.seed(123)
estandarizar(rnorm(100,2,6))
## $med.data
## [1] 2.542435
##
## $m.dv
## [1] -2.352025e-16
##
## $var.x
## [1] 29.99638
##
## $var.dv
## [1] 29.99638
##
## $m.z
## [1] -4.44718e-17
##
## $var.z
## [1] 1
##
## $z
## [1] -0.71304802 -0.35120270 1.60854170 -0.02179795 0.04259548
## [6] 1.77983218 0.40589817 -1.48492941 -0.85149566 -0.58726835
## [11] 1.24195461 0.29513939 0.34000892 0.02221347 -0.70797086
## [16] 1.85854263 0.44636008 -2.25349176 0.66930255 -0.61698896
## [21] -1.26885349 -0.33783464 -1.22304003 -0.89754917 -0.78377819
## [26] -1.94683206 0.81876439 0.06898128 -1.34588242 1.27452758
## [31] 0.36815564 -0.42229479 0.88157948 0.86296437 0.80101058
## [36] 0.65537241 0.50778230 -0.16686565 -0.43422620 -0.51585092
## [41] -0.86009994 -0.32681639 -1.48529653 2.27707482 1.22429519
## [46] -1.32941869 -0.54040552 -0.61026684 0.75541982 -0.19037243
## [51] 0.17847258 -0.13031397 -0.14600575 1.40027842 -0.34637532
## [56] 1.56226982 -1.79571669 0.54141021 0.03664303 0.13752572
## [61] 0.31685861 -0.64934164 -0.46407310 -1.21490140 -1.27319995
## [66] 0.23347834 0.39197814 -0.04097396 0.91131364 2.14685001
## [71] -0.63697081 -2.62876100 1.00275711 -0.87597805 -0.85276181
## [76] 1.02448422 -0.41101270 -1.43635058 0.09957931 -0.25119772
## [81] -0.09272595 0.32303830 -0.50510289 0.60688103 -0.34058618
## [86] 0.26443017 1.10255872 0.37770550 -0.45610238 1.15949090
## [91] 0.98935390 0.50173432 0.16249260 -0.78691881 1.39156928
## [96] -0.75663177 2.29720706 1.57995139 -0.35725306 -1.22349625
### Notese que la desviacion de los datos es igual en los datos
## y en los datos trasladados
N= sort(rnorm(50,4,0.5))
P= sort(rnorm(50,0.5,0.01))
plot(N,P,pch=19, col=2, cex=1.5,ylim=c(0.40,0.6), xlim=c(0,5))
grid(10)
#regracion lineal simple
mod1=lm(P~N)
summary(mod1)
##
## Call:
## lm(formula = P ~ N)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.0034420 -0.0006918 -0.0001520 0.0002590 0.0102823
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4288504 0.0020345 210.79 <2e-16 ***
## N 0.0184706 0.0005212 35.44 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.001805 on 48 degrees of freedom
## Multiple R-squared: 0.9632, Adjusted R-squared: 0.9624
## F-statistic: 1256 on 1 and 48 DF, p-value: < 2.2e-16
abline(a=0.4324012, b=0.01777384, col="black")
abline(h=0.4324012, lty=2, col= "blue")
text(1,0.432, expression(theta), pos=3)
angulo=atan(0.01777384)*180/pi
angulo
## [1] 1.018259
text(1.6,0.432, paste('=',round(angulo,2)),pos=3)
N2= estandarizar(N)$z;N2
## [1] -1.818745675 -1.429286386 -1.378687466 -1.362164698 -1.332455806
## [6] -1.274359642 -1.220877125 -1.219749584 -1.202822125 -1.044268292
## [11] -0.778532245 -0.715588727 -0.705240283 -0.701051654 -0.602227318
## [16] -0.536728829 -0.461427756 -0.402340868 -0.390975754 -0.324912068
## [21] -0.279992511 -0.239208427 -0.170414049 -0.127688018 -0.094651729
## [26] -0.008386507 -0.002215389 0.007286271 0.184560888 0.200476769
## [31] 0.211124580 0.298315209 0.335439103 0.363453249 0.375552705
## [36] 0.494562044 0.516290943 0.561037870 0.713009503 0.781644722
## [41] 0.871156581 0.951971005 0.965988074 1.004562690 1.185542129
## [46] 1.400171993 1.716762417 2.120378565 2.186323503 2.379388117
P2= estandarizar(P)$z;P2
## [1] -1.449671957 -1.395270617 -1.385728059 -1.315058840 -1.218940452
## [6] -1.183848068 -1.168649886 -1.124824437 -1.003432081 -0.971368389
## [11] -0.837984824 -0.698163478 -0.577994948 -0.561333500 -0.534033849
## [16] -0.489448609 -0.444037368 -0.441736445 -0.417258486 -0.390443077
## [21] -0.342868678 -0.295482051 -0.273033477 -0.253479125 -0.169993224
## [26] -0.078277172 -0.001093948 0.028449954 0.049335862 0.059912125
## [31] 0.188660780 0.278654905 0.291816545 0.315148890 0.354636468
## [36] 0.427204421 0.513499503 0.563047008 0.603562908 0.642082578
## [41] 0.768277276 0.784377017 0.804459840 1.007724735 1.089078239
## [46] 1.150528679 1.315157035 2.103608065 2.244576955 3.439655254
plot(N2~P2,pch=19, col=2, cex=1.5)
grid(10)
mod12=lm(P2~N2)
summary(mod12)
##
## Call:
## lm(formula = P2 ~ N2)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.36972 -0.07431 -0.01633 0.02782 1.10447
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.871e-15 2.741e-02 0.00 1
## N2 9.814e-01 2.769e-02 35.44 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1938 on 48 degrees of freedom
## Multiple R-squared: 0.9632, Adjusted R-squared: 0.9624
## F-statistic: 1256 on 1 and 48 DF, p-value: < 2.2e-16
#Calculo de las correlaciones
cor(P,N)
## [1] 0.9814242
cor(P2,N2)
## [1] 0.9814242
Ejercicio sobre fertilización
#Fertilizacion
fert=gl(2,25,50, labels=c("control","quimico"))
datos= data.frame(P,N,P2,N2,fert)
#View(datos)
##prueba T
#h0= Mc=Mq #N
#h0 es falso
#datos originales y estandarizados
#creacion de una funcion para la prueba t
mi.t<- function(x1,x2,conf.lev){
n1=length(x1)
n2=length(x2)
s1=sd(x1)
s2=sd(x2)
tc=(abs((mean(x1)-mean(x2)))/sqrt((((n1-1)*(s1**2)+
(n2-1)*(s2**2))/(n1+n2-2))*((1/n1)+(1+n2))))
p_v=pt(q=tc, df=n1+n2-2, lower.tail = F)
desi= ifelse(p_v < (1-conf.lev)/2,"rechazo Ho"," No rechazo Ho")
return(desi)
}
mi.t(datos$P[1:25],datos$P[26:50],0.95)
## [1] " No rechazo Ho"
mi.t(datos$P2[1:25],datos$P2[26:50],0.95)
## [1] " No rechazo Ho"
#############################################################
(x<-1:5)
## [1] 1 2 3 4 5
x<- 1:5
Ejercicio Ley de gases
ley.gases<- function(v,n,R=0.082,Tp,gra=F){
p=n*R*Tp/v
if(gra) plot3d(Tp, v, p, col=rainbow(100))
return(p)
}
library(rgl)
ley.gases(n=1, Tp=273, v=seq(1,20), gra=T)
## [1] 22.386000 11.193000 7.462000 5.596500 4.477200 3.731000 3.198000
## [8] 2.798250 2.487333 2.238600 2.035091 1.865500 1.722000 1.599000
## [15] 1.492400 1.399125 1.316824 1.243667 1.178211 1.119300
ley.gases(n=1, Tp=rep(seq(273,303,5)), v=rep(seq(1,20),7),gra=T) ## argumentos fijos
## [1] 22.386000 11.398000 7.735333 5.904000 4.805200 4.072667 3.549429
## [8] 2.798250 2.532889 2.320600 2.146909 2.002167 1.879692 1.774714
## [15] 1.492400 1.424750 1.365059 1.312000 1.264526 1.221800 24.846000
## [22] 11.193000 7.598667 5.801500 4.723200 4.004333 3.490857 3.105750
## [29] 2.487333 2.279600 2.109636 1.968000 1.848154 1.745429 1.656400
## [36] 1.399125 1.340941 1.289222 1.242947 1.201300 24.436000 12.423000
## [43] 7.462000 5.699000 4.641200 3.936000 3.432286 3.054500 2.760667
## [50] 2.238600 2.072364 1.933833 1.816615 1.716143 1.629067 1.552875
## [57] 1.316824 1.266444 1.221368 1.180800 24.026000 12.218000 8.282000
## [64] 5.596500 4.559200 3.867667 3.373714 3.003250 2.715111 2.484600
## [71] 2.035091 1.899667 1.785077 1.686857 1.601733 1.527250 1.461529
## [78] 1.243667 1.199789 1.160300 23.616000 12.013000 8.145333 6.211500
## [85] 4.477200 3.799333 3.315143 2.952000 2.669556 2.443600 2.258727
## [92] 1.865500 1.753538 1.657571 1.574400 1.501625 1.437412 1.380333
## [99] 1.178211 1.139800 23.206000 11.808000 8.008667 6.109000 4.969200
## [106] 3.731000 3.256571 2.900750 2.624000 2.402600 2.221455 2.070500
## [113] 1.722000 1.628286 1.547067 1.476000 1.413294 1.357556 1.307684
## [120] 1.119300 22.796000 11.603000 7.872000 6.006500 4.887200 4.141000
## [127] 3.198000 2.849500 2.578444 2.361600 2.184182 2.036333 1.911231
## [134] 1.599000 1.519733 1.450375 1.389176 1.334778 1.286105 1.242300
library(rgl)
##########################
germi= sample(c(0,1), size=100, replace=T);germi
## [1] 1 0 1 1 0 0 1 0 0 1 0 1 0 0 0 0 1 0 0 0 1 1 1 0 1 1 0 0 0 1 1 1 0 1 0
## [36] 1 0 0 0 0 0 1 1 0 1 1 1 0 0 1 1 0 0 0 0 1 1 1 0 0 1 0 0 0 1 0 0 1 0 1
## [71] 0 1 0 1 1 0 1 0 0 1 0 1 0 1 1 0 0 1 0 1 1 0 0 1 0 1 0 0 1 1
pg=(sum(germi)/100)*100
var= replicate(20,sample(c(0,1), size=100, replace=T,prob=c(0.2,0.8)))
colSums(var)
## [1] 73 78 79 87 79 85 82 78 84 75 85 78 72 77 87 81 79 77 82 84
hist(colSums(var))
round(runif(200,0,1)) %>% hist() ###pipes
var2= replicate(200,round(runif(200,0,1)))
pg2=colSums(var2/200)*100
hist(pg2)
Tamaño de muestra
n.muestra.p <- function(N, p, e=0.05, z=1.96){
n=ceiling(N*p*(1-p)/((N-1)*(e/z)^2+p*(1-p)))
}
n=n.muestra.p(N=500, p=0.5);n
## [1] 218