Profesor Oscar Javier Pachecho Perez

Universidad Nacional de Colombia

#“Intervalos de confianza para la media poblaciona, mu” >Suposición 1: La varianza poblacional,sigama^2 es conocida”

set.seed(12)

library(TeachingDemos)

## Warning: package 'TeachingDemos' was built under R version 4.2.2

##install.packages("TeachingDemos")
##?TeachingDemos

POPULATION.MEAN <- 171.66

POPULATION.SD <- 5.60

Repeat at different sample size

##?lapply

#?cat

X <- c(5,10,50,100)

junk <- lapply(X, function(N){
  cat(paste(" sample size ",N))
  ci.examp(mean.sim   = POPULATION.MEAN,
           sd         = POPULATION.SD,
           n          = N,
           reps       = 100,
           conf.level = 0.95,
           method     = "z")

      }
  )

##  sample size  5

##  sample size  10

##  sample size  50

##  sample size  100

##?ci.examp

##install.packages("installr")

library(installr)

## Warning: package 'installr' was built under R version 4.2.2

## 
## Welcome to installr version 0.23.4
## 
## More information is available on the installr project website:
## https://github.com/talgalili/installr/
## 
## Contact: <tal.galili@gmail.com>
## Suggestions and bug-reports can be submitted at: https://github.com/talgalili/installr/issues
## 
##          To suppress this message use:
##          suppressPackageStartupMessages(library(installr))

##updateR() 

#?RGui

#Ejemplo 9.2#

Varianza poblacional es conocida

##install.packages("BSDA")

library(BSDA)

## Warning: package 'BSDA' was built under R version 4.2.2

## Loading required package: lattice

## 
## Attaching package: 'BSDA'

## The following object is masked from 'package:TeachingDemos':
## 
##     z.test

## The following object is masked from 'package:datasets':
## 
##     Orange

n <- 36

xbarra <- 2.6

sigma <- 0.3

Estime mu: La concentración media de zinc en el río

ej1 <- round(zsum.test(mean.x = xbarra,sigma.x = sigma,n.x = n,
                 conf.level = 0.95)$conf.int,2)

ej1

## [1] 2.5 2.7
## attr(,"conf.level")
## [1] 0.95

ej1a <- round(zsum.test(mean.x = xbarra,sigma.x = sigma,n.x = n,
                       conf.level = 0.99)$conf.int,2)
ej1a

## [1] 2.47 2.73
## attr(,"conf.level")
## [1] 0.99

ej1b <- round(zsum.test(mean.x = xbarra,sigma.x = sigma,n.x = n,
                       conf.level = 0.9)$conf.int,2)
ej1b

## [1] 2.52 2.68
## attr(,"conf.level")
## [1] 0.9

Ejemplo 9.3

n ?

Para una confianza del 95

##?qnorm

z.alf.med <- round(qnorm(0.025,mean = 0,sd=1,lower.tail = FALSE),2)

e <- 0.05 # error

sigma <- 0.3 # desviación est conocida (Imposible) 

n <- (z.alf.med*sigma/e)^2

ceiling(n) # Redondea hacia arriba

## [1] 139

cat("el tamaño de la muestra con una confianza del 95%,con sigma conocida 
de 0.3 y un error de 0.03 es",ceiling(n))

## el tamaño de la muestra con una confianza del 95%,con sigma conocida 
## de 0.3 y un error de 0.03 es 139

floor(n) # Redondea hacia arriba

## [1] 138

Sigma es desconocida

Población de donde se extrae la muestra es normal

Si no es normal, con muestra debe ser grande

” Distribución t”

” Cuando la varianza poblacional es desconocida”

pt(2.3,9)-pt(1.5,9)

## [1] 0.06042563

0.9765003-0.9160747

## [1] 0.0604256

qt(0.025,14)

## [1] -2.144787

x <- rt(1000,9)

y <- dt(x,9)

#?plot

plot(x,y,type = "p")


lines(plot(x,y))

P(T>0.8)

grados de libertad, v, es 19

pt(0.8,19,lower.tail = FALSE)

## [1] 0.2167999

Ejemplo 9.5

l <- c(9.8, 10.2, 10.4, 9.8, 10.0, 10.2, 9.6 )

##?t.test
u <- round(t.test(l,conf.level=0.95)$conf.int,2)

u

## [1]  9.74 10.26
## attr(,"conf.level")
## [1] 0.95

paste("Con un 95% de confianza se estima que la media poblacional
      de ácido sulfurico está entre 9.74 y 10.26")

## [1] "Con un 95% de confianza se estima que la media poblacional\n      de ácido sulfurico está entre 9.74 y 10.26"

iter <- stats::rpois(1, lambda = 10)

print an informative message

cat("iteration = ", iter <- iter + 1, "\n")

## iteration =  13

cat(paste(letters, 100* 1:26), fill = TRUE, labels = paste0("{", 1:10, "}:"))

## {1}: a 100 b 200 c 300 d 400 e 500 f 600 g 700 h 800 i 900 j 1000 k 1100 l 1200 
## {2}: m 1300 n 1400 o 1500 p 1600 q 1700 r 1800 s 1900 t 2000 u 2100 v 2200 
## {3}: w 2300 x 2400 y 2500 z 2600

Manual de R

url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'

datos.a <- read.table(file=url, header=T)

View(datos.a)

fix(datos.a)

attach(datos.a)

"Suponga que se quiere obtener un intervalo de confianza bilateral del 
90% para la altura promedio de los hombres de la base de 
datos medidas del cuerpo."

## [1] "Suponga que se quiere obtener un intervalo de confianza bilateral del \n90% para la altura promedio de los hombres de la base de \ndatos medidas del cuerpo."

##install.packages("datos")

library(datos)

## Warning: package 'datos' was built under R version 4.2.2

##?diamantes

fix(diamantes)

hombres <- datos.a[sexo=="Hombre",]

dim(hombres)

## [1] 18  6

fix(hombres)

"Una vez leídos los datos, se analiza la normalidad de la variable altura de los hombres,
a partir de un QQplot y un histograma"

## [1] "Una vez leídos los datos, se analiza la normalidad de la variable altura de los hombres,\na partir de un QQplot y un histograma"

par(mfrow=c(1, 2))

require(car)  # Debe instalar antes el paquete car

## Loading required package: car

## Warning: package 'car' was built under R version 4.2.2

## Loading required package: carData

## Warning: package 'carData' was built under R version 4.2.2

## 
## Attaching package: 'carData'

## The following objects are masked from 'package:BSDA':
## 
##     Vocab, Wool

qqPlot(hombres$altura, pch=19,
       main='QQplot para la altura de hombres',
       xlab='Cuantiles teóricos',
       ylab='Cuantiles muestrales')

## [1] 5 8

hist(hombres$altura, freq=TRUE,
     main='Histograma para la altura de hombres',
     xlab='Altura (cm)',
     ylab='Frecuencia')

"Suponiendo normalidad en las alturas de la base de datos elabore
un intervalo de confianza del 90% para las alturas de los hombres"

## [1] "Suponiendo normalidad en las alturas de la base de datos elabore\nun intervalo de confianza del 90% para las alturas de los hombres"

m1 <- t.test(x=hombres$altura, conf.level=0.90)$conf.int

round(m1,2)

## [1] 176.44 181.72
## attr(,"conf.level")
## [1] 0.9

cat("con una confianza del 90% se puede afirmar que el promedio de las alturas
enlos hombres cae entre",m1)

## con una confianza del 90% se puede afirmar que el promedio de las alturas
## enlos hombres cae entre 176.4384 181.7172

#Intervalo de confianza bilateral para la diferencia de medias poblacionales # de muestras independientes

Ejemplo 2 manual de R

construir un intervalo de confianza del 95%

#para la diferencia de las altura promedio de los hombres y de las mujeres

library(datos)
url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'

datos.a <- read.table(file=url, header=T)

fix(datos.a)

attach(datos.a)

## The following objects are masked from datos.a (pos = 6):
## 
##     altura, biceps, edad, muneca, peso, sexo

hombres <- datos.a[sexo=="Hombre", ]

mujeres <- datos.a[datos.a$sexo=="Mujer", ]

fix(mujeres)

fix(hombres)

Suponer varianzas poblacionales diferentes y

#normalidad en las poblaciones independientes

##?t.test

m2 <- t.test(x=hombres$altura, y=mujeres$altura,
       paired=FALSE, var.equal=FALSE,
       conf.level = 0.95)$conf.int

round(m2,3)

## [1] 10.056 20.033
## attr(,"conf.level")
## [1] 0.95

cat("Con una confianza del 95% se puede afirmar que la diferencia
de promedios poblacionales entre las alturas de hombres y mujeres está entre",m2)

## Con una confianza del 95% se puede afirmar que la diferencia
## de promedios poblacionales entre las alturas de hombres y mujeres está entre 10.05574 20.03315

#Intervalo de confianza bilateral para la #diferencia de medias de muestras pareadas

Ejemplo 9.10 de walpole

Uso Paquete BSDA

n1 <- 50

n2 <- 75

xbarra1 <- 36

xbarra2 <- 42

sigma1 <- 6

sigma2 <- 8

#Librería BSDA

ej4 <- zsum.test(mean.x = xbarra2,mean.y = xbarra1,n.x = n1,
                 n.y = n2,sigma.x = sigma1,sigma.y = sigma2,
                 conf.level = 0.96)$conf.int


ej4

## [1] 3.42393 8.57607
## attr(,"conf.level")
## [1] 0.96

#In de Conf para mB - muA

Intervalo de confianza para muestras pareadas

las diferencias deben ser aprox normales

Ejemplo de manualñ de r

Antes   <- c(81, 87, 86, 82, 90, 86, 96, 73,
             74, 75, 72, 80, 66, 72, 56, 82)
Despues <- c(78, 91, 78, 78, 84, 67, 92, 70,
             58, 62, 70, 58, 66, 60, 65, 73)

Diferencia <- Antes -Despues


require(car)

qqPlot(Diferencia, pch=19, main='QQplot para Diferencias', las=1, 
       xlab='Cuantiles teóricos', ylab='Cuantiles muestrales')

## [1] 15 12

plot(density(Diferencia), main='Densidad para Diferencias', las=1,
     xlab='Diferencia de tiempo', ylab='Densidad')

hist(Diferencia)

ej.mr3 <- t.test(x=Antes, y=Despues, paired=TRUE, conf.level=0.95)$conf.int


round(ej.mr3,3)

## [1]  2.362 11.138
## attr(,"conf.level")
## [1] 0.95

cat("Con una confianmzadel del 95%  se puede afirmar que la diferencia entre el antes y el despue
     cae entre",ej.mr3)

## Con una confianmzadel del 95%  se puede afirmar que la diferencia entre el antes y el despue
##      cae entre 2.362371 11.13763

Intervalo de confianza paera la varianza

Chi-cuadrado

P(X > 10)

##?pchisq

pchisq(10,df=390,lower.tail = FALSE)

## [1] 1

qchisq(0.0470542,19,lower.tail = FALSE)

## [1] 30.38879

set.seed(125)


x <- rchisq(10000,19)

y <- dchisq(x,19)

plot(x,y)

segments(x0=10,y0=0,x1=10,y1=dchisq(10,19))

######## usar paquete stests

if (!require('devtools')) install.packages('devtools')

## Loading required package: devtools

## Warning: package 'devtools' was built under R version 4.2.2

## Loading required package: usethis

## Warning: package 'usethis' was built under R version 4.2.2

devtools::install_github('fhernanb/stests', force=TRUE)

## WARNING: Rtools is required to build R packages, but is not currently installed.
## 
## Please download and install Rtools 4.2 from https://cran.r-project.org/bin/windows/Rtools/ or https://www.r-project.org/nosvn/winutf8/ucrt3/.

## Downloading GitHub repo fhernanb/stests@HEAD

## WARNING: Rtools is required to build R packages, but is not currently installed.
## 
## Please download and install Rtools 4.2 from https://cran.r-project.org/bin/windows/Rtools/ or https://www.r-project.org/nosvn/winutf8/ucrt3/.

## * checking for file 'C:\Users\ACER\AppData\Local\Temp\Rtmp8O1LRW\remotes30883ccc53fe\fhernanb-stests-2760e4f/DESCRIPTION' ... OK
## * preparing 'stests':
## * checking DESCRIPTION meta-information ... OK
## * checking for LF line-endings in source and make files and shell scripts
## * checking for empty or unneeded directories
## Omitted 'LazyData' from DESCRIPTION
## * building 'stests_0.1.0.tar.gz'
## 

## Installing package into 'C:/Users/ACER/AppData/Local/R/win-library/4.2'
## (as 'lib' is unspecified)

require(stests)  # Para cargar el paquete

## Loading required package: stests

## 
## Attaching package: 'stests'

## The following object is masked from 'package:BSDA':
## 
##     z.test

## The following object is masked from 'package:TeachingDemos':
## 
##     z.test

## The following object is masked from 'package:stats':
## 
##     var.test

url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'

datos.a <- read.table(file=url, header=T)

names(datos.a)

## [1] "edad"   "peso"   "altura" "sexo"   "muneca" "biceps"

Ejemplo de manual de r

construir un intervalo de confianza del 98% para la

varianza de la altura de los estudiantes hombres.

fix(hombres)

res <- stests::var.test(x=hombres$altura, conf.level=0.98)

res <- res$conf.int

cat("Con una confianza del 98% se puede afirmar que la varianza poblacional de
    las alturas de los hombres se encuentra entre",res)

## Con una confianza del 98% se puede afirmar que la varianza poblacional de
##     las alturas de los hombres se encuentra entre 21.08468 109.9309

Intervalo de confianza para la proporción

Ejemplo de man ual de r

nc <- 500

xc <- 275

Calcular un intervalo de cfonfianza del 90% para la

prporción de casas con dos televisoreso más

#?prop.test
ej.p <- prop.test(x=xc,n=nc,conf.level = 0.9,correct = TRUE)$conf.int

ej.p <- round(ej.p,3)

ej.p

## [1] 0.512 0.587
## attr(,"conf.level")
## [1] 0.9

cat("Con una confianza del 90% se estima que la proporción de
    casas con mas de dos tv es", ej.p)

## Con una confianza del 90% se estima que la proporción de
##     casas con mas de dos tv es 0.512 0.587

Intervalo de confianza para la diferencia de proporciones

Ejemplo manual r

#Construir un intervalo de confianza del 90% para decidir #si el cambio tuvo efecto positivo o no.

xa <- 75

xd <- 80

na <- 1500

nd <- 2000

ej_dp <- prop.test(x=c(xa,xd),n=c(na,nd),conf.level = 0.9,correct = FALSE)$conf.int

ej_dp <- round(ej_dp,5)

ej_dp

## [1] -0.00173  0.02173
## attr(,"conf.level")
## [1] 0.9

cat("con un nivel de confianza del 90% se estima que la diferencia de proporción de defectuosos 
se encuentra entre",ej_dp)

## con un nivel de confianza del 90% se estima que la diferencia de proporción de defectuosos 
## se encuentra entre -0.00173 0.02173

Otro Ejemolo Walpole

ej.a <- prop.test(x=275, n=500, conf.level=0.90)$conf.int

cat("Con una confianza del 90% se estima que la proporción de casas con más de dos 
    televiores es",ej.a)

## Con una confianza del 90% se estima que la proporción de casas con más de dos 
##     televiores es 0.512231 0.5872162

ejemplo inventado

ej.b <- prop.test(x=c(75,80), n=c(1500,2000),conf.level=0.90)$conf.int

cat("con una confianza del 90% se estima que la diferencia de proporciones de
   los estudiantes que residen afuera entre la unal y la distri está
    ente",ej.b)

## con una confianza del 90% se estima que la diferencia de proporciones de
##    los estudiantes que residen afuera entre la unal y la distri está
##     ente -0.002314573 0.02231457

poblacionales

El estadístco F

Supuesto: las muestras n1 y n2 independientes provienen de poblaciones normales

la distribución maneja grados de libertad en el numerador

grados de libertad en el denominador

gl1= n1-1, gl2 = n2-1

#n1=20, n2=18

set.seed(124)

xf <- rf(n=5000,df1=19,df2=17)

yf <- df(x=xf, df1 = 19,df2 = 17)

plot(xf,yf,main=" Distribución f con gl1=19, gl2=17",col="blue",
     xlab = "F")


area.f <- seq(1.5,6,0.01)

xpf <- c(1.5,area.f,6)

ypf <- c(0,df(area.f,19,17),0)


polygon(xpf,ypf,col = "lightblue")

text(1.8,0.11,"P(F>1.5)=0.2026")

pf(1.5,19,17,lower.tail = FALSE)

## [1] 0.2026496

qf(0.5,19,17,lower.tail = TRUE)

## [1] 1.004227

##3 Ejemplo de manual de r

#Usando la información del ejemplo de diferencia de medias para muestras independientes #se quiere obtener un intervalo de confianza del 95% #para la razón de #las varianzas de las alturas de los estudiantes hombres y mujeres.

url <- 'https://raw.githubusercontent.com/fhernanb/datos/master/medidas_cuerpo'

datos.a <- read.table(file=url, header=T)

hombres <- datos.a[sexo=="Hombre", ]

mujeres <- datos.a[datos.a$sexo=="Mujer", ]

fix(hombres)

fix(mujeres)

inc <- round(var.test(x=hombres$altura,y=mujeres$altura,
         conf.level = 0.95)$conf.int,3)

cat("Con una confianza del 95% se estima que el cociente de varianzas
    poblacionales para las alturas entre hombres y mujeres cae entre",inc)

## Con una confianza del 95% se estima que el cociente de varianzas
##     poblacionales para las alturas entre hombres y mujeres cae entre 0.233 1.663

" se puede suponer que el cociente es igual a 1, lo que indicaría que las
varianzas poblacionales son iguales"

## [1] " se puede suponer que el cociente es igual a 1, lo que indicaría que las\nvarianzas poblacionales son iguales"

Ejercicio

Usar paquuete stests

if (!require('devtools')) install.packages('devtools')

devtools::install_github('fhernanb/stests', force=TRUE)

## WARNING: Rtools is required to build R packages, but is not currently installed.
## 
## Please download and install Rtools 4.2 from https://cran.r-project.org/bin/windows/Rtools/ or https://www.r-project.org/nosvn/winutf8/ucrt3/.

## Downloading GitHub repo fhernanb/stests@HEAD

## WARNING: Rtools is required to build R packages, but is not currently installed.
## 
## Please download and install Rtools 4.2 from https://cran.r-project.org/bin/windows/Rtools/ or https://www.r-project.org/nosvn/winutf8/ucrt3/.

## * checking for file 'C:\Users\ACER\AppData\Local\Temp\Rtmp8O1LRW\remotes308851d45e2c\fhernanb-stests-2760e4f/DESCRIPTION' ... OK
## * preparing 'stests':
## * checking DESCRIPTION meta-information ... OK
## * checking for LF line-endings in source and make files and shell scripts
## * checking for empty or unneeded directories
## Omitted 'LazyData' from DESCRIPTION
## * building 'stests_0.1.0.tar.gz'
## 

## Warning: package 'stests' is in use and will not be installed

library(stests)

stests::var.test  # Para usar la fución del paquete stests

## function (x, y = NULL, alternative = "two.sided", null.value = 1, 
##     conf.level = 0.95) 
## {
##     if (!is.numeric(x)) 
##         stop(paste("The x vector must be numeric", "\n", ""))
##     if (length(x) <= 1) 
##         stop(paste("not enough 'x' observations", "\n", ""))
##     if (!is.null(y)) {
##         if (!is.numeric(y)) 
##             stop(paste("The y vector must be numeric", "\n", 
##                 ""))
##         if (length(y) <= 1) 
##             stop(paste("not enough 'y' observations", "\n", ""))
##     }
##     if (null.value <= 0) 
##         stop(paste("The null value must be positive", "\n", ""))
##     alternative <- match.arg(arg = alternative, choices = c("two.sided", 
##         "greater", "less"))
##     raw.data <- TRUE
##     name_x <- deparse(substitute(x))
##     name_y <- deparse(substitute(y))
##     if (is.null(y)) 
##         res <- var_test_one(var(x), length(x), alternative, conf.level, 
##             null.value, raw.data, name_x)
##     else res <- var_test_two(var(x), length(x), var(y), length(y), 
##         alternative, conf.level, null.value, raw.data, name_x, 
##         name_y)
##     class(res) <- "htest"
##     res
## }
## <bytecode: 0x0000021a410b1938>
## <environment: namespace:stests>

ejercicio

n.w <- 12

n.t <- 10

x.w <- 16

s.w <- 1

x.t <- 11

s.t <- 0.8

ej42 <- var_test(varx = s.w^2,nx=12,vary = s.t^2,ny=n.t,conf.level = 0.98)$conf.int


round(sqrt(ej42),4)

## [1] 0.5493 2.6901
## attr(,"conf.level")
## [1] 0.98