Funciones Estadisticas
Desiderio
12 de Diciembre de 2020
Probabilidad
Funcion que regresa la probabilidad en la cola derecha de Z con media 0 y desviacion 1.
Probabilidad_Z <- function (Valor){
Resultado <- pnorm(Valor, mean = 0, sd = 1, lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de Z dada la probabilidad.
Valor_Z <- function (Probabilidad){
Numero <- qnorm(Probabilidad)
Resultado <- Numero*(-1)
return (Resultado)
}Funcion que regresa la probabilidad de t con x Grados de Libertad.
Probabilidad_t<-function (Numero,Libertad){
Resultado <- pt(Numero,Libertad,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de t dado Alpha y Grados de Libertad.
Valor_t <- function (Probabilidad, Libertad){
Resultado <- qt(Probabilidad,Libertad,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa la probabilidad de Ji Cuadrada con x grados de libertad.
Probabilidad_Ji <- function (Numero, Libertad){
Resultado <- pchisq(Numero,Libertad,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de Ji Cuadrada dado Alpha y Grados de Libertad.
Valor_Ji <- function (Probabilidad, Libertad){
Resultado <- qchisq(Probabilidad,Libertad,lower.tail = FALSE)
return(Resultado)
}Funcion que regresa la probabilidad de F de Fisher dado Numerador y Denominador.
Proba_F <- function (Numero,Numerador,Denominador){
Resultado <- pf(Numero,Numerador,Denominador,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de F de Fisher con 1% de probabilidad dado Numerador y Denominador.
Valor_F_1_Porciento <- function (Numerador,Denominador){
Resultado <- qf(0.01,Numerador,Denominador,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de F de Fisher con 5% de probabilidad dado Numerador y Denominador.
Valor_F_5_Porciento <- function (Numerador,Denominador){
Resultado <- qf(0.05,Numerador,Denominador,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa el valor de F de Fisher con la probabilidad deseada dado Numerador y Denominador.
Valor_F <- function (Probabilidad,Numerador,Denominador){
Resultado <- qf(Probabilidad,Numerador,Denominador,lower.tail = FALSE)
return (Resultado)
}Funcion que regresa la probabilidad Binomial dado Ensayo, Exito y Probabilidad.
Probabilidad_Binomial <- function (Ensayo,Exito,Probabilidad){
R <- dbinom(Exito,Ensayo,Probabilidad,log=FALSE)
return (R)
}Funcion que regresa la probabilidad de Poisson dado Lambda y K.
Probabilidad_Poisson <- function (K,Lambda){
R <- dpois(K,Lambda)
return (R)
}Funcion que suma todos los elementos de un vector.
Suma <- function (x){
z<-0
for (i in x){
z<-z+i
}
return (z)
}Funcion que calcula la suma cuadrada de los elementos de un vector.
Suma_Cuadrada <- function(x){
z<-0
for (i in x){
a<-i^2
z<-z+a
}
return (z)
}Funcion que calcula la suma del producto de los numeros de dos vectores.
Suma_xy <- function (x,y){
Indicador <- 0
Limite <- length(x)
for (i in 1:Limite){
A<-x[i]
B<-y[i]
Indicador<-Indicador+(A*B)
}
return (Indicador)
}Funcion que calcula la desviacion estandar poblacional de un vector.
Desviacion_Estandar_Poblacional <- function (Vector){
Promedio <- mean(Vector)
n <- length(Vector)
p <- 0
for (i in Vector){
a <- i
Resta <- a-Promedio
Resta_C <- Resta^(2)
p <- p + Resta_C
}
p <- p^(1/2)
n <- n^(1/2)
return (p/n)
}Funcion que calcula la varianza poblacional de un vector.
Varianza_Poblacional <- function (Vector){
Resultado <- Desviacion_Estandar_Poblacional(Vector)^2
return (Resultado)
}Funcion que calcula el valor absoluto de un numero.
Valor_Absoluto <- function (Numero){
if (Numero<0){
Numero <- Numero*(-1)
return (Numero)
}else{
return (Numero)
}
}SE DEFINEN LOS VECTORES QUE SERAN USADOS EN LOS EJEMPLOS.
x <- c(4,6,8,10,12,14,16)
y <- c(1,3,6,8,14,16,20)Regresión Lineal
Funcion que calcula Sxx.
Sxx <- function (Vector){
Longitud <- length(Vector)
Suma_Cuadrados <- Suma_Cuadrada(Vector)
Suma_X <- Suma(Vector)
Suma_X_Cuadrada <- Suma_X^2
Division <- Suma_X_Cuadrada/Longitud
return (Suma_Cuadrados-Division)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Sxx
Sxx(x)## [1] 112Funcion que calcula Syy
Syy <- function (Vector){
Longitud <- length(Vector)
Suma_Cuadrados <- Suma_Cuadrada(Vector)
Suma_Y <- Suma(Vector)
Suma_Y_Cuadrada <- Suma_Y^2
Division <- Suma_Y_Cuadrada/Longitud
return (Suma_Cuadrados-Division)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Syy
Syy(y)## [1] 301.4286Funcion que calcula Sxy
Sxy <- function (X,Y){
Longitud <- length(X)
S_XY <- Suma_xy(X,Y)
S_X <- Suma(X)
S_Y <- Suma(Y)
Producto <- S_X*S_Y
Division <- Producto/Longitud
return (S_XY-Division)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Sxy
Sxy(x,y)## [1] 182Funcion que calcula el coeficiente de regresion Betta
Betta <- function (x,y){
sxy <- Sxy(x,y)
sxx <- Sxx(x)
return (sxy/sxx)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Betta
Betta(x,y)## [1] 1.625Funcion que calcula el coeficiente de regresion Alpha
Alpha <- function (x,y){
B<-Betta(x,y)
Media_X<-mean(x)
Media_Y<-mean(y)
Producto<-B*Media_X
return (Media_Y-Producto)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Alpha
Alpha(x,y)## [1] -6.535714Funcion que calcula los coeficientes de regresion
Coeficientes_Regresion <- function (x,y){
A<-Alpha(x,y)
B<-Betta(x,y)
R<-c(A,B)
return (R)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Coeficientes_Regresion
Coeficientes_Regresion(x,y)## [1] -6.535714 1.625000#NOTEMOS QUE EL PRIMER VALOR REGRESADO CORRESPONDE A Alpha Y EL SEGUNDO A BettaFuncion que calcula la ecuacion linea de minimos cuadrados ajustada a los datos
Ecuacion_Linea_Minimos_Cuadrados <- function (X,Y){
Alpha <- Coeficientes_Regresion(X,Y)[1]
Betta <- Coeficientes_Regresion(X,Y)[2]
Mensaje1<-"y="
Mensaje2<-"x"
Mensaje3<-"+"
Resultado <- c(Mensaje1,Alpha,Mensaje3,Betta,Mensaje2)
return(Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Ecuacion Linea de Minimos Cuadrados
Ecuacion_Linea_Minimos_Cuadrados(x,y)## [1] "y=" "-6.53571428571429" "+"
## [4] "1.625" "x"Funcion que calcula el valor de y dado x, los elementos de x, los elementos de y
Estimacion_de_y_Min_Cuadrados <- function (Valor_X,X,Y){
A<-Alpha(X,Y)
B<-Betta(X,Y)
Producto<-B*Valor_X
return (A+Producto)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Estimacion_de_y_Min_Cuadrados
Estimacion_de_y_Min_Cuadrados(7,x,y)## [1] 4.839286Estadística
NOTA: Varianza_Poblacional solo debe tomar dos valores: TRUE o FALSE, el primero es en caso de conocer la varianza poblacional y el segundo es en caso de que la desconozca.
Funcion que calcula el intervalo de confianza para una media conociendo la media de la muestra, desviacion, tamaño, probabilidad deseada y varianza poblacional.
Intervalo_Confianza_Una_Media_Datos_Conocidos <- function (Media, Desviacion, n, Probabilidad, Varianza_Poblacional) {
Mensaje <- " < M < "
Raiz_n<-n^(1/2)
Division<-Desviacion/Raiz_n
if (n>=30 || Varianza_Poblacional == TRUE){
ValorZ <- Valor_Z(Probabilidad)
Producto<-Division*ValorZ
S<-Media+Producto
R<-Media-Producto
Resultado <- c(R,Mensaje,S)
return (Resultado)
}else if(n<=29 && Varianza_Poblacional == FALSE){
t<-Valor_t(Probabilidad,n-1)
Producto <- Division*t
S<-Media+Producto
R<-Media-Producto
Resultado<-c(R,Mensaje,S)
return (Resultado)
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para una media con datos conocidos
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.1
n <- 20
Media <- 64.3
Desviacion <- 15
Probabilidad <- 0.025
Varianza_Poblacional <- TRUE
Intervalo_Confianza_Una_Media_Datos_Conocidos(Media,Desviacion,n,Probabilidad,Varianza_Poblacional)## [1] "57.7260809456756" " < M < " "70.8739190543244"Funcion que calcula el intervalo de confianza para una media desconociendo media de la muestra,desviacion, dado un vector con los numeros, probabilidad y varianza poblacional.
Intervalo_Confianza_Una_Media_Datos_Desconocidos <- function (Vector,Probabilidad,Varianza_Poblacional){
Mensaje <- "< M <"
n<-length(Vector)
Raiz_n<-n^(1/2)
Media<-mean(Vector)
Desviacion<-sd(Vector)
if (n>=30 || Varianza_Poblacional == TRUE){
Desviacion_Poblacional <- Desviacion_Estandar_Poblacional(Vector)
Division <- Desviacion_Poblacional/Raiz_n
ValorZ <- Valor_Z(Probabilidad)
Producto <- ValorZ*Division
S<-Media+Producto
R<-Media-Producto
Resultado <- c(R,Mensaje,S)
return (Resultado)
}else if (n<=29 && Varianza_Poblacional == FALSE){
Desviacion_Muestral <- sd(Vector)
Division <- Desviacion_Muestral/Raiz_n
Valort <- Valor_t(Probabilidad,n-1)
Producto <- Valort*Division
S <- Media+Producto
R <- Media-Producto
Resultado <- c(R,Mensaje,S)
return (Resultado)
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para una media con datos desconocidos
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.3
Numeros <- c(17,13,18,19,17,21,29,22,16,28,21,15,26,23,24,20,8,17,17,21,32,18,25,22,16,10,20,22,19,14,30,22,12,24,28,11)
Probabilidad <- 0.025
Varianza_Poblacional <- FALSE
Intervalo_Confianza_Una_Media_Datos_Desconocidos(Numeros,Probabilidad,Varianza_Poblacional)## [1] "18.0697977075648" "< M <" "21.7635356257685"Funcion que calcula el intervalo de confianza para diferencia de medias conociendo la media de ambos grupos, desviacion, probabilidad y Varianza poblacional.
Intervalo_Confianza_Diferencia_Media_Datos_Conocidos <- function (Media1,Desviacion1,n1,Media2,Desviacion2,n2,Varianza_Poblacional1,Varianza_Poblacional2,Probabilidad){
Mensaje <- "< M1 - M2 <"
Diferencia_Medias<-Media1-Media2
if (n1>=30 && n2>=30 && Varianza_Poblacional1 == TRUE && Varianza_Poblacional2 == TRUE){
ValorZ <- Valor_Z(Probabilidad)
Varianza1<-Desviacion1^2
Varianza2<-Desviacion2^2
Division1<-Varianza1/n1
Division2<-Varianza2/n2
Suma_Division<-Division1+Division2
Raiz_Suma_Division<-Suma_Division^(1/2)
Producto<-Raiz_Suma_Division*ValorZ
S<-Diferencia_Medias+Producto
R<-Diferencia_Medias-Producto
Resultado<-c(R,Mensaje,S)
if (R<0 && S>0){
print("El 0 se encuentra en el intervalo de confianza.")
print("Las medias son iguales")
return (Resultado)
}else if (R<0 && S<0){
print("El intervalo de confianza es negativo.")
print("Por lo tanto M2>M1")
return (Resultado)
}else if (R>0 && S>0){
print("El intervalo de confianza es positivo.")
print("Por lo tanto M1>M2")
return (Resultado)
}
}else if (n1<=29 && n2<=29 && Varianza_Poblacional1 == FALSE && Varianza_Poblacional2 == FALSE){
Valort <- Valor_t(Probabilidad,n1+n2-2)
Varianza1<-Desviacion1^2
Varianza2<-Desviacion2^2
#CALCULAMOS Sp, INICIO
X<-(n1-1)*Varianza1
Y<-(n2-1)*Varianza2
Suma_Numerador<-X+Y
Suma_Denonimador<-n1+n2-2
Division_Sp<-Suma_Numerador/Suma_Denonimador
Sp<-Division_Sp^(1/2)
#CALCULAMOS Sp, FIN
A<-1/n1
B<-1/n2
C<-A+B
R<-C^(1/2)
Prod<-Valort*Sp*R
Suma<-Diferencia_Medias+Prod
Resta<-Diferencia_Medias-Prod
Resultado<-c(Resta,Mensaje,Suma)
if (Resta<0 && Suma>0){
print("El 0 se encuentra en el intervalo de confianza.")
print("Las medias son iguales.")
return (Resultado)
}else if(Resta<0 && Suma<0){
print("El intervalo de confianza es negativo.")
print("Por lo tanto M2>M1")
return (Resultado)
}else if(Resta>0 && Suma>0){
print("El intervalo de confianza es positivo.")
print("Por lo tanto M1>M2")
return (Resultado)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Diferencia de Medias con Datos Conocidos
#LOS DATOS SE OBTUVIERON DEL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO: https://www.youtube.com/watch?v=vcV_-KSzy0M&t=45s
Sindicalizada_Salario_Medio <- 20.75
Sindicalizada_Desviacion_Estandar <- 2.25
Sindicalizada_Tamaño <- 40
No_Sindicalizada_Salario_Medio <- 19.80
No_Sindicalizada_Desviacion_Estandar <- 1.90
No_Sindicalizada_Tamaño <- 45
Intervalo_Confianza_Diferencia_Media_Datos_Conocidos(Sindicalizada_Salario_Medio,Sindicalizada_Desviacion_Estandar,Sindicalizada_Tamaño,No_Sindicalizada_Salario_Medio,No_Sindicalizada_Desviacion_Estandar,No_Sindicalizada_Tamaño,TRUE,TRUE,0.05)## [1] "El intervalo de confianza es positivo."
## [1] "Por lo tanto M1>M2"## [1] "0.202026035469634" "< M1 - M2 <" "1.69797396453036"Funcion que calcula el intervalo de confianza para diferencia de medias con datos desconocidos con base a probabilidad y Varianza Poblacional.
Intervalo_Confianza_Diferencia_Media_Datos_Desconocidos <- function (Vector1,Vector2,Varianza_Poblacional1,Varianza_Poblacional2,Probabilidad){
Mensaje <- "< M1 - M2 <"
n1 <- length(Vector1)
n2 <- length(Vector2)
Media1<-mean(Vector1)
Media2<-mean(Vector2)
Diferencia_Medias<-Media1-Media2
if (n1>=30 && n2>=30 && Varianza_Poblacional1 == TRUE && Varianza_Poblacional2 == TRUE){
ValorZ<-Valor_Z(Probabilidad)
Varianza1<-Varianza_Poblacional(Vector1)
Varianza2<-Varianza_Poblacional(Vector2)
Division1<-Varianza1/n1
Division2<-Varianza2/n2
Suma_Division<-Division1+Division2
Raiz_Suma_Division<-Suma_Division^(1/2)
Producto<-Raiz_Suma_Division*ValorZ
Suma<-Diferencia_Medias+Producto
Resta<-Diferencia_Medias-Producto
Resultado<-c(Resta,Mensaje,Suma)
if (Resta>0 && Suma>0){
print("El intervalo de confianza es positivo")
print("Por lo tanto M1>M2")
return (Resultado)
}else if (Resta<0 && Suma<0){
print("El intervalo de confianza es negativo")
print("Por lo tanto M2>M1")
return (Resultado)
}else if (Resta<0 && Suma >0){
print("El 0 esta en el intervalo de confianza")
print("Por lo tanto M1 = M2")
return (Resultado)
}
}else if (n1<=29 && n2<=29 && Varianza_Poblacional1 == FALSE && Varianza_Poblacional2 == FALSE){
t<-Valor_t(Probabilidad,n1+n2-2)
Varianza1<-var(Vector1)
Varianza2<-var(Vector2)
#CALCULAMOS Sp, INICIO
X<-(n1-1)*Varianza1
Y<-(n2-1)*Varianza2
Suma_Numerador<-X+Y
Suma_Denonimador<-n1+n2-2
Division_Sp<-Suma_Numerador/Suma_Denonimador
Sp<-Division_Sp^(1/2)
#CALCULAMOS Sp, FIN
A<-1/n1
B<-1/n2
C<-A+B
R<-C^(1/2)
Prod<-t*Sp*R
Suma<-Diferencia_Medias+Prod
Resta<-Diferencia_Medias-Prod
Resultado<-c(Resta,Mensaje,Suma)
if (Resta>0 && Suma>0){
print("El intervalo de confianza es positivo")
print("Por lo tanto M1>M2")
return (Resultado)
}else if (Resta<0 && Suma<0){
print("El intervalo de confianza es negativo")
print("Por lo tanto M2>M1")
return (Resultado)
}else if (Resta<0 && Suma >0){
print("El 0 esta en el intervalo de confianza")
print("Por lo tanto M1 = M2")
return (Resultado)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Diferencia de Medias con Datos Desconocidos
#LOS DATOS SE OBTUVIERON DEL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO https://www.youtube.com/watch?v=YIImQIRDGMg&t=198s
Urbano <- c(3,7,25,10,15,6,12,25,15,7)
Rural <- c(48,44,40,38,33,21,20,12,1,18)
Intervalo_Confianza_Diferencia_Media_Datos_Desconocidos(Urbano,Rural,FALSE,FALSE,0.05)## [1] "El intervalo de confianza es negativo"
## [1] "Por lo tanto M2>M1"## [1] "-24.4006223702352" "< M1 - M2 <" "-5.59937762976485"Funcion que calcula el intervalo de confianza para proporciones dado el numero de total, n, y numero muestral, x.
Intervalo_Confianza_Proporciones <- function (N,X,Probabilidad){
Mensaje <- "< Tetta <"
if (X>=30){
Tetta<-X/N
Z<-Valor_Z(Probabilidad)
Numerador<-Tetta*(1-Tetta)
Division<-Numerador/N
Raiz_Cuadrada<-Division^(1/2)
Producto<-Z*Raiz_Cuadrada
Suma<-Tetta+Producto
Resta<-Tetta-Producto
Resultado<-c(Resta,Mensaje,Suma)
return (Resultado)
}else{
print("El tamaño de la muestra es menor a 30.")
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Proporciones
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.7
n<-400
x<-136
probabilidad<-0.025
Intervalo_Confianza_Proporciones(n,x,probabilidad)## [1] "0.293577393455242" "< Tetta <" "0.386422606544758"Funcion que calcula el intervalo de confianza para una diferencia de proporciones dado el numero total y muestral de ambos grupos.
Intervalo_Confianza_Diferencia_Proporciones <- function (N1,X1,N2,X2,Probabilidad){
Mensaje<-"< Tetta1 - Tetta2 <"
if (X1>=30 && X2>=30){
Tetta1<-X1/N1
Tetta2<-X2/N2
Diferencia_Tetta <- Tetta1-Tetta2
Z<-Valor_Z(Probabilidad)
Parte1<-(Tetta1*(1-Tetta1))/N1
Parte2<-(Tetta2*(1-Tetta2))/N2
Suma_0<-Parte1+Parte2
Raiz<-Suma_0^(1/2)
Producto<-Z*Raiz
Suma<-Diferencia_Tetta+Producto
Resta<-Diferencia_Tetta-Producto
Resultado<-c(Resta,Mensaje,Suma)
if (Resta<0 && Suma<0){
print("El intervalo de confianza es negativo")
print("Por lo tanto Tetta2>Tetta1")
return(Resultado)
}else if (Resta>0 && Suma>0){
print("El intervalo de confianza es positivo")
print("Por lo tanto Tetta1>Tetta2")
return (Resultado)
}else if (Resta<0 && Suma>0){
print("El 0 esta en el intervalo de confianza")
print("Por lo tanto Tetta1 = Tetta2")
return (Resultado)
}
}else{
print("El tamaño de la muestra es menor a 30.")
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Diferencia de Proporciones
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.9
Total_Hombres <- 200
Muestra_Hombres <- 132
Total_Mujeres <- 159
Muestra_Mujeres <- 90
Probabilidad <- 0.005
Intervalo_Confianza_Diferencia_Proporciones(Total_Hombres,Muestra_Hombres,Total_Mujeres,Muestra_Mujeres,Probabilidad)## [1] "El 0 esta en el intervalo de confianza"
## [1] "Por lo tanto Tetta1 = Tetta2"## [1] "-0.0390587954610169" "< Tetta1 - Tetta2 <" "0.226983323762904"Funcion que calcula el intervalo de confianza para varianza dado el tamaño, varianza y probabilidad.
Intervalo_Confianza_Varianza_Datos_Conocidos <- function (n,varianza,probabilidad){
Mensaje <- "< Sigma^2 <"
Ji_Izquierda<-Valor_Ji(probabilidad,n-1)
Ji_Derecha<-Valor_Ji(1-probabilidad,n-1)
Numerador<-(n-1)*varianza
Izquierda<-Numerador/Ji_Izquierda
Derecha<-Numerador/Ji_Derecha
Resultado<-c(Izquierda,Mensaje,Derecha)
return (Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Varianza con Datos Conocidos
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.10
n <- 16
Desviacion <- 2.2
Varianza <- Desviacion^2
Probabilidad <- 0.005
Intervalo_Confianza_Varianza_Datos_Conocidos(n,Varianza,Probabilidad)## [1] "2.21332551771247" "< Sigma^2 <" "15.7794679924421"Funcion que calcula el intervalo de confianza para varianza con datos desconocidos.
Intervalo_Confianza_Varianza_Datos_Desconocidos <- function (Vector,Probabilidad){
Mensaje <- "< Sigma^2 <"
n<-length(Vector)
Varianza<-var(Vector)
Ji_Iz<-Valor_Ji(Probabilidad,n-1)
Ji_De<-Valor_Ji(1-Probabilidad,n-1)
Numerador<-(n-1)*Varianza
Iz<-Numerador/Ji_Iz
De<-Numerador/Ji_De
Resultado<-c(Iz,Mensaje,De)
return (Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Varianza con Datos Desconocidos
#LOS DATOS SE OBTUVIERON DEL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO https://www.youtube.com/watch?v=Pw66ZHxzAF4
Datos<-c(6.0,6.4,7.0,5.8,6.0,5.8,5.9,6.7,6.1,6.5,6.3,5.8)
Probabilidad <- 0.025
Intervalo_Confianza_Varianza_Datos_Desconocidos(Datos,Probabilidad)## [1] "0.0770603499359095" "< Sigma^2 <" "0.442682943162402"Funcion que calcula el intervalo de confianza de dos varianzas dado tamaño y varianza de ambos grupos.
Intervalo_Confianza_Cociente_Varianzas_Datos_Conocidos <- function (n1,varianza1,n2,varianza2,probabilidad){
Mensaje <- "< (Sigma1^2)/(Sigma2^2) <"
C_V<-varianza1/varianza2
F_Iz<-Valor_F(probabilidad,n1-1,n2-1)
F_De<-Valor_F(probabilidad,n2-1,n1-1)
Division<-1/F_Iz
Izquierda<-C_V*Division
Derecha<-C_V*F_De
Resultado<-c(Izquierda,Mensaje,Derecha)
if (Izquierda<1 && Derecha>1){
print("El intervalo de confianza contiene al 1")
print("Por lo tanto Sigma1^2 == Sigma2^2")
return(Resultado)
}
return (Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Cociente de Varianzas con Datos Conocidos.
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 11.11
n1 <- 10
n2 <- 8
Desviacion1 <- 0.5
Desviacion2 <- 0.7
Varianza1 <- Desviacion1^2
Varianza2 <- Desviacion2^2
Probabilidad <- 0.01
Intervalo_Confianza_Cociente_Varianzas_Datos_Conocidos(n1,Varianza1,n2,Varianza2,Probabilidad)## [1] "El intervalo de confianza contiene al 1"
## [1] "Por lo tanto Sigma1^2 == Sigma2^2"## [1] "0.0759373236345369" "< (Sigma1^2)/(Sigma2^2) <"
## [3] "2.86370687621237"Funcion que calcula el intervalo de confianza de dos varianzas con datos desconocidos.
Intervalo_Confianza_Cociente_Varianzas_Datos_Desconocidos <- function (Vector1,Vector2,Probabilidad){
Mensaje <- "< (Sigma1^2)/(Sigma2^2) <"
n1<-length(Vector1)
n2<-length(Vector2)
Varianza1<-var(Vector1)
Varianza2<-var(Vector2)
C_V<-Varianza1/Varianza2
F_Iz<-Valor_F(Probabilidad,n1-1,n2-1)
F_De<-Valor_F(Probabilidad,n2-1,n1-1)
Division<-1/F_Iz
Izquierda<-C_V*Division
Derecha<-C_V*F_De
Resultado<-c(Izquierda,Mensaje,Derecha)
if (Izquierda<1 && Derecha>1){
print("El intervalo de confianza contiene al 1")
print("Por lo tanto Sigma1^2 == Simga2^2")
return (Resultado)
}
return (Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Intervalo de Confianza para Cociente de Varianza con Datos Conocidos
#LOS DATOS SE OBTUVIERON EN EL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO https://www.youtube.com/watch?v=SzXfVJbPh-0
Gasolina <- c(34,36,39,31,33,26,45,34,39,38,37)
Gas <- c(33,41,39,32,29,28,33,34,25,28,36,33,35)
Probabilidad <- 0.05
Intervalo_Confianza_Cociente_Varianzas_Datos_Desconocidos(Gasolina,Gas,Probabilidad)## [1] "El intervalo de confianza contiene al 1"
## [1] "Por lo tanto Sigma1^2 == Simga2^2"## [1] "0.436251043034655" "< (Sigma1^2)/(Sigma2^2) <"
## [3] "3.49897398512978"NOTA: PARA LAS PRUEBA DE HIPOTESIS SE TIENEN TRES OPCIONES AL ARGUMENTO ‘Hipotesis’.
Opcion “1”: Algo >= Algo Opcion “2”: Algo <= Algo Opcion “3”: Algo != Algo
ADEMAS, INTRODUCIR LA PROBABILIDAD EN SU FORMA DECIMAL, ES DECIR, SI REQUIERES EL 90% ENTONCES INTRODUCES 0.1, PERO SI LA FORMULA INDICA Alpha/2 ENTONCES INTRODUCES 0.05, LAS FUNCIONES NO CALCULAN EL “Alpha/2” NI CONVIERTEN LOS PORCENTAJES A DECIMAL.
Funcion que calcula la prueba de hipotesis para una media.
Prueba_Hipotesis_Una_Media <- function (n,Media,Desviacion,Mo,Probabilidad,Varianza_Poblacional,Hipotesis){
M1 <- "Rechazamos Ho"
M2 <- "No rechazamos Ho"
Estadistico_Prueba <- (Media-Mo)/((Desviacion)/(n^(1/2)))
if ((n>=30 && Varianza_Poblacional == TRUE) || (n>=30 && Varianza_Poblacional == FALSE) && (Hipotesis == 1)){
Z <- Valor_Z(Probabilidad)
if (Estadistico_Prueba>=Z){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de Z es")
print(Z)
return (M2)
}
}else if (((n>=30 && Varianza_Poblacional == TRUE) || (n>=30 && Varianza_Poblacional == FALSE)) && (Hipotesis == 2)){
Z <- Valor_Z(Probabilidad)
ZN <- Z*(-1)
if (Estadistico_Prueba<=ZN){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de Z Negativo es")
print(ZN)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de Z Negativo es")
print(ZN)
return (M2)
}
}else if (((n>=30 && Varianza_Poblacional == TRUE) || (n>=30 && Varianza_Poblacional == FALSE)) && (Hipotesis == 3)){
Z <- Valor_Z(Probabilidad)
if (Valor_Absoluto(Estadistico_Prueba)>=Z){
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de Z es")
print(Z)
return (M2)
}
}else if ((n<=29 && Varianza_Poblacional == FALSE) && Hipotesis == 1){
t <- Valor_t(Probabilidad,n-1)
if (Estadistico_Prueba>=t){
return (M1)
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de t es")
print(t)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de t es")
print(t)
return (M2)
}
}else if ((n<=29 && Varianza_Poblacional == FALSE) && Hipotesis == 2){
t <- Valor_t(Probabilidad,n-1)
tN <- t*(-1)
if (Estadistico_Prueba<=tN){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de t Negativo es")
print(tN)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba)
print("El valor de t Negativo es")
print(tN)
return (M2)
}
}else if ((n<=29 && Varianza_Poblacional == FALSE) && Hipotesis == 3){
t <- Valor_t(Probabilidad,n-1)
if (Valor_Absoluto(Estadistico_Prueba)>=t){
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de t es")
print(t)
return (M1)
}else{
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de t es")
print(t)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis Una Media
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 13.2
n <- 100
Media <- 21819
Desviacion <- 1295
Mo <- 22000
Probabilidad <- 0.05
Varianza_Poblacional <- FALSE
Hipotesis <- 2
Prueba_Hipotesis_Una_Media(n,Media,Desviacion,Mo,Probabilidad,Varianza_Poblacional,Hipotesis)## [1] "El Estadistico de Prueba es"
## [1] -1.397683
## [1] "El valor de Z Negativo es"
## [1] -1.644854## [1] "No rechazamos Ho"Funcion que calcula la prueba de hipotesis para dos medias.
#NOTA: SI SE DESEA SABER SI UNA MEDIA ES MAYOR A OTRA, DARLE VALOR 0 A Mo
Prueba_Hipotesis_Dos_Media <- function (n1,Media1,Varianza1,n2,Media2,Varianza2,Mo,Probabilidad,Varianza_Poblacional,Hipotesis){
M1 <- "Rechazamos Ho"
M2 <- "No rechazamos Ho"
#CALCULAMOS EL ESTADISTICO DE PRUEBA PARA LA Z, INICIO
Numerador <- Media1-Media2-Mo
A<-Varianza1/n1
B<-Varianza2/n2
Suma_Denominador<-A+B
Denominador<-Suma_Denominador^(1/2)
Estadistico_Prueba_Z<-Numerador/Denominador
#CALCULAMOS EL ESTADISTICO DE PRUEBA PARA LA Z, FIN
#CALCULAMOS EL ESTADISTICO DE PRUEBA PARA LA t, INICIO
Num<-((n1-1)*Varianza1)+((n2-1)*Varianza2)
Den<-n1+n2-2
Div<-Num/Den
Sp<-Div^(1/2)
A<-(1/n1)+(1/n2)
Raiz2<-A^(1/2)
D<-Sp*Raiz2
Estadistico_Prueba_t<-Numerador*(1/D)
#CALCULAMOS EL ESTADISTICO DE PRUEBA PARA LA t, FIN
if ( ((n1>=30 && n2>=30) && (Varianza_Poblacional==TRUE)) || (n1>=30 && n2>=30) && Hipotesis == 1 ){
Z<-Valor_Z(Probabilidad)
if (Estadistico_Prueba_Z>=Z){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_Z)
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_Z)
print("El valor de Z es")
print(Z)
return (M2)
}
}else if (((n1>=30 && n2>=30) && (Varianza_Poblacional==TRUE)) || (n1>=30 && n2>=30) && Hipotesis == 2){
Z<-Valor_Z(Probabilidad)
ZN<-Z*(-1)
if (Estadistico_Prueba_Z<=ZN){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_Z)
print("El valor de Z Negativo es")
print(ZN)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_Z)
print("El valor de Z Negativo es")
print(ZN)
return (M2)
}
}else if (((n1>=30 && n2>=30) && (Varianza_Poblacional==TRUE)) || (n1>=30 && n2>=30) && Hipotesis == 3){
Z<-Valor_Z(Probabilidad)
if (Valor_Absoluto(Estadistico_Prueba_Z)>=Z){
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba_Z))
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba_Z))
print("El valor de Z es")
print(Z)
return (M2)
}
}else if (n1<=29 && n2<=29 && Varianza_Poblacional == FALSE && Hipotesis == 1){
t<-Valor_t(Probabilidad,n1+n2-2)
if (Estadistico_Prueba_t>=t){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_t)
print("El valor de t es")
print(t)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_t)
print("El valor de t es")
print(t)
return (M2)
}
}else if (n1<=29 && n2<=29 && Varianza_Poblacional == FALSE && Hipotesis == 2){
t<-Valor_t(Probabilidad,n1+n2-2)
tN<-t*(-1)
if (Estadistico_Prueba_t<=tN){
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_t)
print("El valor de t Negativo es")
print(tN)
return (M1)
}else{
print("El Estadistico de Prueba es")
print(Estadistico_Prueba_t)
print("El valor de t Negativo es")
print(tN)
return (M2)
}
}else if (n1<=29 && n2<=29 && Varianza_Poblacional == FALSE && Hipotesis == 3){
t<-Valor_t(Probabilidad,n1+n2-2)
if (Valor_Absoluto(Estadistico_Prueba_t)>=t){
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba_t))
print("El valor de t es")
print(t)
return (M1)
}else{
print("El valor absoluto del Estadistico de Prueba es")
print(Valor_Absoluto(Estadistico_Prueba_t))
print("El valor de t es")
print(t)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis para Dos Medias
#LOS DATOS SE OBTUVIERON EN EL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO https://www.youtube.com/watch?v=8K7QbatdrIA
Mujeres <- c(12,11,14,13,13,14,13,12,14,12)
Hombres <- c(13,10,11,12,13,12,10,12)
Mujeres_n <- length(Mujeres)
Hombres_n <- length(Hombres)
Mujeres_Media <- mean(Mujeres)
Hombres_Media <- mean(Hombres)
Mujeres_Varianza <- var(Mujeres)
Hombres_Varianza <- var(Hombres)
Mo <- 0
Probabilidad <- 0.05
Varianza_Poblacional <- FALSE
Hipotesis <- 1 #QUEREMOS PROBAR QUE LA MEDIA DE LAS MUJERES ES MAYOR AL DE LOS HOMBRES
Prueba_Hipotesis_Dos_Media(Mujeres_n,Mujeres_Media,Mujeres_Varianza,Hombres_n,Hombres_Media,Hombres_Varianza,Mo,Probabilidad,Varianza_Poblacional,Hipotesis)## [1] "El Estadistico de Prueba es"
## [1] 2.245266
## [1] "El valor de t es"
## [1] 1.745884## [1] "Rechazamos Ho"Funcion que calcula la prueba de hipotesis para una varianza.
#VARIANZA0 HACE REFERENCIA AL DATO QUE ESTAMOS COMPARANDO.
Prueba_Hipotesis_Una_Varianza <- function (n,Varianza,Probabilidad,Varianza0,Hipotesis){
M1<-"Rechazamos Ho"
M2<-"No rechazamos Ho"
Estadistico_Prueba <- ((n-1)*Varianza)/(Varianza0)
if (Hipotesis == 1){
Ji<-Valor_Ji(Probabilidad,n-1)
if (Estadistico_Prueba>=Ji){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M2)
}
}else if (Hipotesis == 2){
Ji<-Valor_Ji(1-Probabilidad,n-1)
if (Estadistico_Prueba<=Ji){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M2)
}
}else if (Hipotesis == 3){
Ji<-Valor_Ji(Probabilidad,n-1)
if (Estadistico_Prueba>=Ji){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis Una Varianza
#LOS DATOS SE OBTUVIERON EN EL EJERCICIO REALIZADO EN EL SIGUIENTE VIDEO https://www.youtube.com/watch?v=qrtBrueW_sE
n <- 20
Varianza <- 22.63
Probabilidad <- 0.01
Varianza0 <- 18
Hipotesis <- 1 #QUEREMOS PROBAR QUE LA VARIANZA ES MAYOR A 18
Prueba_Hipotesis_Una_Varianza(n,Varianza,Probabilidad,Varianza0,Hipotesis)## [1] "El estadistico de prueba es"
## [1] 23.88722
## [1] "El valor de Ji es"
## [1] 36.19087## [1] "No rechazamos Ho"Funcion que calcula la prueba de hipotesis para dos varianzas.
Prueba_Hipotesis_Dos_Varianza <- function (n1,Varianza1,n2,Varianza2,Probabilidad,Hipotesis){
M1 <- "Rechazamos Ho"
M2 <- "No rechazamos Ho"
Estadistico_Prueba1 <- Varianza1/Varianza2
Estadistico_Prueba2 <- Varianza2/Varianza1
if (Hipotesis == 1){
f<-Valor_F(Probabilidad,n1-1,n2-1)
if (Estadistico_Prueba1>=f){
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
print("El valor de f es")
print(f)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
print("El valor de f es")
print(f)
return (M2)
}
}else if (Hipotesis == 2){
f<-Valor_F(Probabilidad,n2-1,n1-1)
if (Estadistico_Prueba2>=f){
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
print("El valor de f es")
print(f)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
print("El valor de f es")
print(f)
return (M2)
}
}else if (Hipotesis == 3){
if (Varianza1>Varianza2){
f<-Valor_F(Probabilidad, n1-1,n2-1)
if (Estadistico_Prueba1>=f){
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
print("El valor de f es")
print(f)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
print("El valor de f es")
print(f)
return (M2)
}
}else if (Varianza2>Varianza1){
f<-Valor_F(Probabilidad,n2-1,n1-1)
if (Estadistico_Prueba2>=f){
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
print("El valor de f es")
print(f)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
print("El valor de f es")
print(f)
return (M2)
}
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis para Dos Varianzas
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 13.7
n1 <- 13
Varianza1 <- 19.2
n2 <- 16
Varianza2 <- 3.5
Probabilidad <- 0.02
Hipotesis <- 3 #QUEREMOS PROBAR QUE LAS VARIANZAS SON DISTINTAS
Prueba_Hipotesis_Dos_Varianza(n1,Varianza1,n2,Varianza2,Probabilidad,Hipotesis)## [1] "El estadistico de prueba es"
## [1] 5.485714
## [1] "El valor de f es"
## [1] 3.127944## [1] "Rechazamos Ho"Funcion que calcula la prueba de hipotesis para proporciones (muestra chica)
FUNCION EXTRA 1. Funcion que calcula el estadistico de prueba 1.
E_P_1 <- function (n,Tetta0,Probabilidad){
Contador<-0
for (i in n:1){
A<-Probabilidad_Binomial(n,i,Tetta0)
Contador<-Contador+A
if (Contador>Probabilidad){
break
}
}
return (i+1)
}FUNCION EXTRA 2. Funcion que calcula el estadistico de prueba 2.
E_P_2 <- function (n,Tetta0,Probabilidad){
Contador<-0
for (i in 0:10000){
A<-Probabilidad_Binomial(n,i,Tetta0)
Contador<-Contador+A
if (Contador>Probabilidad){
break
}
}
return (i-1)
}Prueba_Hipotesis_Proporciones_Chico <- function (n,x,Tetta0,Probabilidad,Hipotesis){
M1<-"Rechazamos Ho"
M2<-"No rechazamos Ho"
if (Hipotesis == 1){
#CALCULAMOS EL ESTADISTICO DE PRUEBA, INICIO
Estadistico_Prueba1<-E_P_1(n,Tetta0,Probabilidad)
#CALCULAMOS EL ESTADISTICO DE PRUEBA, FIN
if (x>=Estadistico_Prueba1){
print("El numero de Exitos es")
print(x)
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
return (M1)
}else{
print("El numero de Exitos es")
print(x)
print("El estadistico de prueba es")
print(Estadistico_Prueba1)
return (M2)
}
}else if (Hipotesis == 2){
#CALCULAMOS EL ESTADISTICO DE PRUEBA, INICIO
Estadistico_Prueba2<-E_P_2(n,Tetta0,Probabilidad)
#CALCULAMOS EL ESTADISTICO DE PRUEBA, FIN
if (x<=Estadistico_Prueba2){
print("El numero de Exitos es")
print(x)
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
return (M1)
}else{
print("El numero de Exitos es")
print(x)
print("El estadistico de prueba es")
print(Estadistico_Prueba2)
return (M2)
}
}else if (Hipotesis == 3){
#CALCULAMOS EL ESTADISTICO DE PRUEBA 1
E_P1 <- E_P_1(n,Tetta0,Probabilidad)
#CALCULAMOS EL ESTADISTICO DE PRUEBA 2
E_P2 <- E_P_2(n,Tetta0,Probabilidad)
if ((x>=E_P1) || (x<=E_P2)){
print("El numero de exitos es")
print(x)
print("El estadistico de prueba 1 es")
print(Estadistico_Prueba1)
print("El estadistico de prueba 2 es")
print(Estadistico_Prueba2)
return (M1)
}else{
print("El numero de exitos es")
print(x)
print("El estadistico de prueba 1 es")
print(Estadistico_Prueba1)
print("El estadistico de prueba 2 es")
print(Estadistico_Prueba2)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis para Proporciones con muestra chica
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJERCICIO 13.54
Repeticiones <- 19
Exitos <- 1
Tetta0 <- 0.3
Probabilidad <- 0.05
Hipotesis <- 2 #QUEREMOS PROBAR QUE Tetta ES MENOR A 0.30
Prueba_Hipotesis_Proporciones_Chico(Repeticiones,Exitos,Tetta0,Probabilidad,Hipotesis)## [1] "El numero de Exitos es"
## [1] 1
## [1] "El estadistico de prueba es"
## [1] 2## [1] "Rechazamos Ho"Funcion que calcula la prueba de hipotesis para proporciones (Muestra grande)
Prueba_Hipotesis_Proporciones_Grande <- function (n,x,Tetta0,Probabilidad,Hipotesis){
M1<-"Rechazamos Ho"
M2<-"No rechazamos Ho"
Numerador<-x-(n*Tetta0)
A<-(n*Tetta0)*(1-Tetta0)
Denominador<-A^(1/2)
Estadistico_Prueba<-Numerador/Denominador
if (Hipotesis == 1){
Z<-Valor_Z(Probabilidad)
if (Estadistico_Prueba>=Z){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Z es")
print(Z)
return (M2)
}
}else if (Hipotesis == 2){
Z<-Valor_Z(Probabilidad)
ZN<-Z*(-1)
if (Estadistico_Prueba<=ZN){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Z Negativo es")
print(ZN)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Z Negativo es")
print(ZN)
return (M2)
}
}else if (Hipotesis == 3){
Z<-Valor_Z(Probabilidad)
if (Valor_Absoluto(Estadistico_Prueba)>=Z){
print("El valor absoluto del estadistico de prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de Z es")
print(Z)
return (M1)
}else{
print("El valor absoluto del estadistico de prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de Z es")
print(Z)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis para proporciones con muestra grande.
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJERCICIO 13.62
n <- 200
Exitos <- 174
Tetta0 <- 0.9
Probabilidad <- 0.05
Hipotesis <- 2 #QUEREMOS PROBAR QUE Tetta ES MENOR A 0.9
Prueba_Hipotesis_Proporciones_Grande(n,Exitos,Tetta0,Probabilidad,Hipotesis)## [1] "El estadistico de prueba es"
## [1] -1.414214
## [1] "El valor de Z Negativo es"
## [1] -1.644854## [1] "No rechazamos Ho"Funcion que calcula la prueba de hipotesis para K proporciones.
NOTA: Si desconoces el valor de Tetta0 escribe FALSE. En los argumentos tienes que ingresar los valores correspondientes al Exito y Fracaso acomodados, es decir, correspondientes a cada valor.
Prueba_Hipotesis_K_Proporciones <- function (Exito,Fracaso,Probabilidad,Tetta0){
M1="Rechazamos Ho"
M2="No rechazamos Ho"
Matriz <- cbind(Exito,Fracaso)
colnames(Matriz)<-c("Exito","Fracaso")
K<-length(Exito)
Repeticiones<-c()
for (i in 1:K){
A<-Exito[i]
B<-Fracaso[i]
Suma<-A+B
Repeticiones<-c(Repeticiones,Suma)
}
Matriz<-cbind(Matriz,Repeticiones)
if (Tetta0 != FALSE){
Ji<-Valor_Ji(Probabilidad,K)
Exitos_Esperados<-c()
for (j in 1:K){
A1<-Repeticiones[j]
Producto<-A1*Tetta0
Exitos_Esperados<-c(Exitos_Esperados,Producto)
}
Fracasos_Esperados<-c()
for (ñ in 1:K){
A2<-Repeticiones[ñ]
Produc<-A2*(1-Tetta0)
Fracasos_Esperados<-c(Fracasos_Esperados,Produc)
}
Matriz<-cbind(Matriz,Exitos_Esperados,Fracasos_Esperados)
print("La matriz con los exitos observados, fracasos observados, repeticiones, exitos esperados y fracasos esperados es la siguiente")
print(Matriz)
Pivote1<-0
for (u in 1:K){
Elemento1<-Exito[u]
Elemento2<-Exitos_Esperados[u]
Resta_U<-Elemento1-Elemento2
Resta_U_C<-Resta_U^2
Division_U<-Resta_U_C/Elemento2
Pivote1<-Pivote1+Division_U
}
Pivote2<-0
for (y in 1:K){
Elemento3<-Fracaso[y]
Elemento4<-Fracasos_Esperados[y]
Resta_Y<-Elemento3-Elemento4
Resta_Y_C<-Resta_Y^2
Division_Y<-Resta_Y_C/Elemento4
Pivote2<-Pivote2+Division_Y
}
Estadistico_Prueba<-Pivote1+Pivote2
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de Ji es")
print(Ji)
if (Estadistico_Prueba>=Ji){
return (M1)
}else{
return (M2)
}
}else if (Tetta0 == FALSE){
Ji<-Valor_Ji(Probabilidad,K-1)
G<-Suma(Exito)
D<-Suma(Repeticiones)
Tetta0<-G/D
Exitos_Esperados<-c()
for (j in 1:K){
A1<-Repeticiones[j]
Producto<-A1*Tetta0
Exitos_Esperados<-c(Exitos_Esperados,Producto)
}
Fracasos_Esperados<-c()
for (ñ in 1:K){
A2<-Repeticiones[ñ]
Produc<-A2*(1-Tetta0)
Fracasos_Esperados<-c(Fracasos_Esperados,Produc)
}
Matriz<-cbind(Matriz,Exitos_Esperados,Fracasos_Esperados)
print("La matriz con los exitos observados, fracasos observados, repeticiones, exitos esperados y fracasos esperados es la siguiente")
print(Matriz)
Pivote1<-0
for (u in 1:K){
Elemento1<-Exito[u]
Elemento2<-Exitos_Esperados[u]
Resta_U<-Elemento1-Elemento2
Resta_U_C<-Resta_U^2
Division_U<-Resta_U_C/Elemento2
Pivote1<-Pivote1+Division_U
}
Pivote2<-0
for (y in 1:K){
Elemento3<-Fracaso[y]
Elemento4<-Fracasos_Esperados[y]
Resta_Y<-Elemento3-Elemento4
Resta_Y_C<-Resta_Y^2
Division_Y<-Resta_Y_C/Elemento4
Pivote2<-Pivote2+Division_Y
}
Estadistico_Prueba<-Pivote1+Pivote2
print("El estadistico de prueba es")
print(as.character(Estadistico_Prueba))
print("El valor de Ji es")
print(Ji)
if (Estadistico_Prueba>=Ji){
return (M1)
}else{
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Prueba de Hipotesis para K Proporciones.
#EN UN ESTUDIO SOBRE UN TALLER, SE OBTUVO UN CONJUNTO DE DATOS PARA DETERMINAR SI LA PROPORCION
#DE ARTICULOS DEFECTUOSO PRODUCIDOS POR LOS OBREROS ERA O NO LA MISMA PARA TURNOS DIURNO,
#VESPERTINO Y NOCTURNO. SE OBTUVIERON LOS SIGUIENTES RESULTADOS.
Defectuosos <- c(45,55,70)
No_Defectuosos <- c(905,890,870)
Matriz_Ejemplo_K <- cbind(Defectuosos,No_Defectuosos)
rownames(Matriz_Ejemplo_K)<-c("Diurno","Vespertino","Nocturno")
Matriz_Ejemplo_K## Defectuosos No_Defectuosos
## Diurno 45 905
## Vespertino 55 890
## Nocturno 70 870#USAR UN NIVEL DE SIGNIFICANCIA DE 0.05 PARA DETERMINAR SI LA PROPORCION DE ARTICULOS DEFECTUOSOS ES LA MISMA
Exito <- Matriz_Ejemplo_K[,1]
Fracasos <- Matriz_Ejemplo_K[,2]
Probabilidad <- 0.05
Tetta0 <- FALSE #NO CONOCEMOS EL VALOR DE Tetta0
Prueba_Hipotesis_K_Proporciones(Exito,Fracasos,Probabilidad,Tetta0)## [1] "La matriz con los exitos observados, fracasos observados, repeticiones, exitos esperados y fracasos esperados es la siguiente"
## Exito Fracaso Repeticiones Exitos_Esperados Fracasos_Esperados
## Diurno 45 905 950 56.96649 893.0335
## Vespertino 55 890 945 56.66667 888.3333
## Nocturno 70 870 940 56.36684 883.6332
## [1] "El estadistico de prueba es"
## [1] "6.23392004768477"
## [1] "El valor de Ji es"
## [1] 5.991465## [1] "Rechazamos Ho"Funcion que calcula la prueba de hipotesis para tablas RxC.
NOTA: Esta funcion recibe como parametro la matriz que se desea analizar y la probabilidad deseada.
Tablas_RxC <- function (Matriz,Probabilidad){
M1<-"Rechazamos Ho"
M2<-"No rechazamos Ho"
Numero_Filas<-nrow(Matriz)
Numero_Columnas<-ncol(Matriz)
f_i <- c()
f_j <- c()
for (i in 1:Numero_Filas){
Fila_i <- Matriz[i,]
Suma_Fila_i <- Suma(Fila_i)
f_i<-c(f_i,Suma_Fila_i)
}
for (j in 1:Numero_Columnas){
Columna_j <- Matriz[,j]
Suma_Columna_j <- Suma(Columna_j)
f_j <-c(f_j,Suma_Columna_j)
}
f<-Suma(f_j)
Frecuencias_Esperadas <- c()
for (u in 1:Numero_Filas){
A<-f_i[u]
for (p in 1:Numero_Columnas){
B<-f_j[p]
Resultado<-(A*B)/f
Frecuencias_Esperadas<-c(Frecuencias_Esperadas,Resultado)
}
}
Vector1<-c()
for (m in 1:Numero_Columnas){
Objeto<-Matriz[m,]
Vector1<-c(Vector1,Objeto)
}
Limite<-length(Vector1)
Pivote<-0
for (v in 1:Limite){
Numero1<-Vector1[v]
Numero2<-Frecuencias_Esperadas[v]
Resta<-Numero1-Numero2
Resta_C<-Resta^2
Division <- Resta_C/Numero2
Pivote <- Pivote + Division
}
Estadistico_Prueba<-Pivote
Ji<-Valor_Ji(Probabilidad,(Numero_Columnas-1)*(Numero_Filas-1))
if (Estadistico_Prueba>=Ji){
print("El estadistico de prueba es")
print(as.character(Estadistico_Prueba))
print("El valor de Ji es")
print(Ji)
return (M1)
}else{
print("El estadistico de prueba es")
print(as.character(Estadistico_Prueba))
print("El valor de Ji es")
print(Ji)
return (M2)
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Tablas RxC
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJERCICIO 13.11
Baja <- c(63,58,14)
Promedio <- c(42,61,47)
Altas <- c(15,31,29)
Matriz_Calificaciones <- cbind(Baja,Promedio,Altas)
rownames(Matriz_Calificaciones)<-c("Bajo","Promedio","Alto")
Matriz_Calificaciones## Baja Promedio Altas
## Bajo 63 42 15
## Promedio 58 61 31
## Alto 14 47 29Probabilidad <- 0.01
Tablas_RxC(Matriz_Calificaciones,Probabilidad)## [1] "El estadistico de prueba es"
## [1] "32.1398518518519"
## [1] "El valor de Ji es"
## [1] 13.2767## [1] "Rechazamos Ho"Funcion que calcula la prueba de Bondad y Ajuste.
ANTES DE PODER USAR ESTA FUNCION, ES NECESARIO INSTALAR LA PAQUETERIA ‘MASS’ y ‘survival’ y ‘fitdistrplus’.
Tome en cuenta que la funcion solo ejecuta funciones dentro de las paqueterias y hace mas ameno al programador el uso de estas, sin embargo, hay un paso que no es posible programar de manera general y es el nombre del archivo CSV y columna que se desea analizar, por lo tanto, se dejan las variables necesarias para su trabajo y facil acceso. Se describen a continuacion los valores a ingresar en el parametro de la funcion para determinar que distribucion se desea analizar.
Binomial <-> 1 Poisson <-> 2 Geometrica <-> 3 Hipergeometrica <-> 4 Binomial Negativa <-> 5 Uniforme <-> 6 Normal <-> 7 t Student <-> 8 Chi-Cuadrado <-> 9 Exponencial <-> 10 Gamma <-> 11 Weibull <-> 12
library(MASS)
library(survival)
library(fitdistrplus)
# Elimina los # una vez que hayas terminado de ingresar los datos necesarios.
# Archivo<-read.csv("Escribe aqui el nombre del archivo csv que deseas analizar, recuerda agregar la terminacion .csv")
# Dato<-Archivo$#Escoge el nombre de la columna que deseas analizar los datos
Prueba_Bondad_Ajuste <- function (Dato,Distribucion) {
if (Distribucion == 1){
Ajuste <- fitdist(Dato,"binom")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 2){
Ajuste <- fitdist(Dato,"pois")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 3){
Ajuste <- fitdist(Dato,"geom")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 4){
Ajuste <- fitdist(Dato,"hyper")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 5){
Ajuste <- fitdist(Dato,"nbinom")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 6){
Ajuste <- fitdist(Dato,"unif")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 7){
Ajuste <- fitdist(Dato,"norm")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 8){
Ajuste <- fitdist(Dato,"t")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 9){
Ajuste <- fitdist(Dato,"chisq")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 10){
Ajuste <- fitdist(Dato,"exp")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 11){
Ajuste <- fitdist(Dato,"gamma")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}else if (Distribucion == 12){
Ajuste <- fitdist(Dato,"weibull")
Ajuste <- gofstat(Ajuste)
return(Ajuste$kstest)
}
}
#EN ESTE APARTADO SOLO MODIFICA EL VALOR EN BASE A LA DISTRIBUCION QUE DESEAS ANALIZAR.
#SE IMPRIMIRA LAS GRAFICAS TEORICAS Y SE IMPRIMIRA UN MENSAJE QUE INDICA SI RECHAZAMOS O NO LA HIPOTESIS NULA.
#Prueba_Bondad_Ajuste(Ingresa el numero correspondiente a la funcion analizar,Dato)SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION DE PRUEBA DE BONDAD Y AJUSTE.
Quiero analizar si mi funcion 1 se distribuye de manera Uniforme y mi funcion 2 se distribuye de manera exponencial.
Funcion1 <- read.csv("Datos_Funcion1.csv")
Funcion2 <- read.csv("Datos_Funcion2.csv")Sospecho que mi hipotesis es cierta pues la grafica de la funcion 1 y funcion 2 respectivamente son las siguientes.
library(ggplot2)
Grafo_Funcion1 <- ggplot(Funcion1)
Grafo_Funcion1 <- Grafo_Funcion1 + aes( x = Primos, y = Tiempo)
Grafo_Funcion1 <- Grafo_Funcion1 + geom_point()
Grafo_Funcion1 <- Grafo_Funcion1 + ggtitle("Grafica Funcion 1")
Grafo_Funcion1
Grafo_Funcion2 <- ggplot(Funcion2)
Grafo_Funcion2 <- Grafo_Funcion2 + aes( x = Primos, y = Tiempo)
Grafo_Funcion2 <- Grafo_Funcion2 + geom_point()
Grafo_Funcion2 <- Grafo_Funcion2 + ggtitle("Grafica Funcion 2")
Grafo_Funcion2
Realizo la prueba de bondad y ajuste.
Prueba de Bondad y Ajuste para la funcion 1.
Prueba_Bondad_Ajuste(Funcion1$Tiempo,6)## 1-mle-unif
## "rejected"Prueba de Bondad y Ajuste para la funcion 2.
Prueba_Bondad_Ajuste(Funcion2$Tiempo,10)## 1-mle-exp
## "rejected"Funcion que realiza el analisis de regresion lineal.
Analisis_Regresion_Lineal <- function (Vector_X,Vector_Y,Bo,Probabilidad,Hipotesis) {
M1<-"Rechazamos Ho"
M2<-"No rechazamos Ho"
n<-length(Vector_X)
Sxx <- Sxx(Vector_X)
Syy <- Syy(Vector_Y)
Sxy <- Sxy(Vector_X,Vector_Y)
B <- Betta(Vector_X,Vector_Y)
#CALCULAMOS EL ESTADISTICO DE PRUEBA, INICIO
Sigma_Cuadrada<-(Syy-(B*Sxy))/n
Sigma<-Sigma_Cuadrada^(1/2)
A_t<-(B-Bo)/Sigma
Bt_0<-((n-1)*Sxx)/n
B_t<-Bt_0^(1/2)
Estadistico_Prueba<-A_t+B_t
#CALCULAMOS EL ESTADISTICO DE PRUEBA, FIN
if (Hipotesis == 1){
t <- Valor_t(Probabilidad,n-2)
if (Estadistico_Prueba>=t){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de t es")
print(t)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de t es")
print(t)
return (M2)
}
}else if (Hipotesis == 2){
t <- Valor_t(Probabilidad,n-2)
tN<-t*(-1)
if (Estadistico_Prueba<=tN){
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de t Negativo es")
print(tN)
return (M1)
}else{
print("El estadistico de prueba es")
print(Estadistico_Prueba)
print("El valor de t Negativo es")
print(tN)
return (M2)
}
}else if (Hipotesis == 3){
t <- Valor_t(Probabilidad,n-2)
if (Valor_Absoluto(Estadistico_Prueba)>=t){
print("El valor absoluto del estadistico de prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de t es")
print(t)
return (M1)
}else{
print("El valor absoluto del estadistico de prueba es")
print(Valor_Absoluto(Estadistico_Prueba))
print("El valor de t es")
print(t)
return (M2)
}
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Analisis de Regresion Lineal
#LOS DATOS SE OBTUVIERON DE EJERCICIOS PERTENECIENTES AL LIBRO 'ESTADISTICA MATEMATICA CON APLICACIONES', SEXTA EDICION, JOHN. FREUND, IRWIN MILLER, MARYLESS MYLER. EJEMPLO 14.5
x<-c(4,9,10,14,4,7,12,22,1,17)
y<-c(31,57,65,73,37,44,60,91,24,84)
Bo<-3
Probabilidad <- 0.01
Hipotesis <- 1 #QUEREMOS VERIFICAR QUE B>3
Analisis_Regresion_Lineal(x,y,Bo,Probabilidad,Hipotesis)## [1] "El estadistico de prueba es"
## [1] 18.48642
## [1] "El valor de t es"
## [1] 2.896459## [1] "Rechazamos Ho"Funcion que calcula el Estimador de Maxima Verosimilitud p De la correlacion entre dos variables aleatorias X y Y.
Estimador_Max_Verosimilitud_p <- function (Vector_X,Vector_Y){
Sxx <- Sxx(Vector_X)
Syy <- Syy(Vector_Y)
Sxy <- Sxy(Vector_X,Vector_Y)
D_1 <- Sxx*Syy
Denominador <- D_1^(1/2)
print("El estimador de maxima verosimilitud es")
return (Sxy/Denominador)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Estimador Maxima Verosimilitud
#SE DIERON DIVERSAS DOSIS DE UNA SUSTANCIA VENENOSA A GRUPOS DE 25 RATONES
#Y SE OBSERVARON LOS SIGUIENTES RESULTADOS
X<-c(4,6,8,10,12,14,16)
Y<-c(1,3,6,8,14,16,20)
Matriz_1 <- cbind(X,Y)
colnames(Matriz_1)<-c("Dosis","Numero de muertes")
Matriz_1## Dosis Numero de muertes
## [1,] 4 1
## [2,] 6 3
## [3,] 8 6
## [4,] 10 8
## [5,] 12 14
## [6,] 14 16
## [7,] 16 20#ESTIMAMOS EL ESTIMADOR DE MAXIMA VEROSIMILITUD, p, DE LA CORRELACION ENTRE DOS VARIABLES ALEATORIAS X Y Y.
r<-Estimador_Max_Verosimilitud_p(X,Y)## [1] "El estimador de maxima verosimilitud es"r## [1] 0.9905358Funcion que determina el porcentaje de variacion total de Y.
Porcentaje_Variacion_Total_Y <- function (r){
r_cuadrado<-r^(2)
Porcentaje_Variacion<-100*r_cuadrado
if (Porcentaje_Variacion == 1){
Mensaje<-"Todos los puntos caen sobre la recta pues el valor del porcentaje de variacion es 1"
return (Mensaje)
}else if (Porcentaje_Variacion>=80 && Porcentaje_Variacion<=99.99){
Valor<-as.character(Porcentaje_Variacion)
Mensaje<-"El ajuste lineal es idoneo pues el porcentaje de variacion es "
Resultado<-c(Mensaje,Valor)
return (Resultado)
}else if (Porcentaje_Variacion>=50 && Porcentaje_Variacion<=79.99){
Valor<-as.character(Porcentaje_Variacion)
Mensaje<-"El ajuste lineal no es un representacion real de los datos, existe una perdida de datos, se recomienda usar otro ajuste pues el porcentaje de variacion es"
Resultado<-c(Mensaje,Valor)
return(Resultado)
}else if (Porcentaje_Variacion>=0 && Porcentaje_Variacion<=49.99){
Valor<-as.character(Porcentaje_Variacion)
Mensaje<-"El ajuste lineal no es idoneo, usa otro ajuste pues el porcentaje de variacion es "
Resultado<- c(Mensaje,Valor)
return (Resultado)
}
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Porcentaje Variacion Total de Y
#UTILIZAREMOS EL EJEMPLO ANTERIOR PARA REALIZAR EL PRESENTE
Porcentaje_Variacion_Total_Y(r)## [1] "El ajuste lineal es idoneo pues el porcentaje de variacion es "
## [2] "98.1161137440758"Funcion que determina el estimador de Maxima Verosimilitud p y el porcentaje de variacion total de Y.
Estimador_Maxima_Verosimilitud_P_Porcentaje_Variacion_Y <- function (Vector_X,Vector_Y){
r <- Estimador_Max_Verosimilitud_p(Vector_X,Vector_Y)
Mensaje<-Porcentaje_Variacion_Total_Y(r)
Resultado <- c(r,Mensaje)
return (Resultado)
}SE ANEXA UN EJEMPLO PARA EL USO DE LA FUNCION Estimador Maxima Verosimilitud p y Porcentaje de Variacion Y
#UTILIZAREMOS EL EJEMPLO DE LA FUNCION DE ESTIMADOR DE MAXIMA VEROSIMILITUD P
Estimador_Maxima_Verosimilitud_P_Porcentaje_Variacion_Y(X,Y)## [1] "El estimador de maxima verosimilitud es"## [1] "0.990535783018846"
## [2] "El ajuste lineal es idoneo pues el porcentaje de variacion es "
## [3] "98.1161137440758"Analisis de Regresión
El siguiente apartado contiene funciones que diseñe en mi clase de Analisis de Regresion. Estas funciones trabajan con matrices pues son parte de una regresion lineal multivariada.
La teoria detras de estas funciones y aplicaciones la puede encontrar en:
Kutner, M. Nachtsheim, C. Neter, J. Li, W. (2004). Applied Linear Statistical Models. (5). McGraw-Hill Irwin.
Funcion para diseñar la matriz X
#EN CASO DE TENER MAS VARIABLES Xi AGREGAR A LOS ARGUMENTOS DE LA FUNCION Y AL VECTOR M DENTRO DE LA FUNCION
M_X <- function (Vector_X1,Vector_X2,Vector_X3){
Vector <- c()
for (i in 1:length(Vector_X1)){
Vector <- c(Vector,1)
}
M <- cbind(Vector,Vector_X1,Vector_X2,Vector_X3)
return (M)
}Funcion para calcular la matriz H
H <- function (M_X){
T_M_X <- t(M_X)
A1 <- T_M_X%*%M_X
A <- solve(A1)
return (M_X%*%A%*%T_M_X)
}Funcion para calcular los estimadores de las B_i
#EN CASO DE TENER MAS VARIABLES Xi AGREGAR MAS bi EN EL NOMBRE DE LAS FILAS
Matriz_b <- function(M_X,M_Y){
T_M_X <- t(M_X)
A1 <- T_M_X%*%M_X
A <- solve(A1)
B <- T_M_X%*%M_Y
Resultado <- A%*%B
row.names(Resultado) <- c("b0","b1","b2","b3")
colnames(Resultado) <- c("Coeficientes de Regresion")
return (Resultado)
}Funcion que calcula los valores ajustados
Matriz_Valores_Ajustados <- function (M_X,M_Y){
Resultado <- H(M_X)%*%M_Y
colnames(Resultado) <- c("Valores ajustados")
return (Resultado)
}Funcion que calcula los residuos
Matriz_residuos <- function (M_X,M_Y){
Matriz_H <- H(M_X)
Matriz_Identidad <- diag(1, nrow = nrow(Matriz_H), ncol = ncol(Matriz_H))
A <- Matriz_Identidad-Matriz_H
Resultado <- A%*%M_Y
colnames(Resultado) <- c("Residuos")
return (Resultado)
}Funcion para diseñar la matriz de valores Xh
#LA CANTIDAD DE ARGUMENTOS EN LA FUNCION DEPENDE LA CANTIDAD DE VARIABLES DADAS EN LA MATRIZ X
Matriz_Xh <- function (X_h1,X_h2,X_h3){
Vector <- c(1,X_h1,X_h2,X_h3)
Vector2 <- cbind(Vector)
colnames(Vector2) <- c("Valores Xh")
return (Vector2)
}Funcion que calcula la respuesta media estimada
Matriz_Gorro_Yh <- function (M_X,M_Y,M_Xh){
b <- Matriz_b(M_X,M_Y)
a <- t(M_Xh)
resultado <- a%*%b
colnames(resultado) <- c("Respuesta media estimada")
return(resultado)
}Funciones para calcular el ANOVA
SSTO <- function (M_Y){
n <- length(M_Y)
Matriz_Identidad <- diag(1, nrow = n, ncol = n)
Matriz_1 <- matrix(1,nrow = n, ncol = n)
Transpuesta <- t(M_Y)
A <- (Matriz_Identidad-(1/n)*Matriz_1)
return (Transpuesta%*%A%*%M_Y)
}
SSE <- function (M_X,M_Y){
n <- length(M_Y)
Matriz_H <- H(M_X)
Matriz_Identidad <- diag(1, nrow = nrow(Matriz_H), ncol = ncol(Matriz_H))
Transpuesta <- t(M_Y)
A <- Matriz_Identidad - Matriz_H
return (Transpuesta%*%A%*%M_Y)
}
SSR <- function (M_X,M_Y){
n <- length(M_Y)
Matriz_H <- H(M_X)
Matriz_1 <- matrix(1, nrow = n, ncol= n)
Transpuesta <- t(M_Y)
A <- (Matriz_H - (1/n)*Matriz_1)
return (Transpuesta%*%A%*%M_Y)
}Funcion que calcula el ANOVA
ANOVA <- function (M_X,M_Y,Probabilidad){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
ssr <- SSR(M_X,M_Y)
sse <- SSE(M_X,M_Y)
ssto <- SSTO(M_Y)
msr <- ssr/p_1
mse <- sse/n_p
Estadistico_Prueba <- msr/mse
Critico <- Valor_F(Probabilidad,p_1,n_p)
FUENTE_VARIACION <- c("Regresion","Error","Total")
SS <- c(ssr,sse,ssto)
GL <- c(p_1,n_p,n-1)
MS <- c(msr,mse,"-")
ESTADISTICO_PRUEBA <- c(Estadistico_Prueba,"-","-")
VALOR_CRITICO <- c(Critico,"-","-")
M <- cbind(FUENTE_VARIACION,SS,GL,MS,ESTADISTICO_PRUEBA,VALOR_CRITICO)
return (M)
}Funcion que calcula la matriz de covarianza del residuo
Matriz_Varianza_Covarianza_de_Matriz_Residuo <- function (M_X,M_Y){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
mse <- SSE(M_X,M_Y)/n_p
Matriz_H <- H(M_X)
Matriz_Identidad <- diag(1, nrow = nrow(Matriz_H), ncol = ncol(Matriz_H))
A <- Matriz_Identidad-Matriz_H
return (mse[1]*A)
}Funcion que determina el coeficiente de determinacion
COEFICIENTE_DETERMINACION <- function (M_X,M_Y){
ssr <- SSR(M_X,M_Y)
ssto <- SSTO(M_Y)
resultado <- ssr/ssto
colnames(resultado) <- c("Coeficiente de determinacion")
row.names(resultado) <- c("[1]")
return (resultado)
}Funcion que determina la contribucion de una variable a un modelo
Contribucion <- function (Matriz_X_Original,Matriz_X_Nueva,y){
E_Full <- ANOVA(Matriz_X_Original,y,0.05)[2,2]
E_Red <- ANOVA(Matriz_X_Nueva,y,0.05)[2,2]
return (as.numeric(E_Full) - as.numeric(E_Red))
}Funcion que calcula el S_b
# EL NUMERO DE COLUMNAS Y FILAS DEPENDE DEL NUMERO DE COLUMNAS Y FILAS QUE HAYAN EN LA MATRIZ X, SI TENEMOS, POR EJEMPLO, UNA MATRIZ 3X3 ENTONCES DAMOS 4 VALORES A LAS FILAS Y COLUMNAS.
S_b <- function (M_X,M_Y){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
sse <- SSE(M_X,M_Y)
mse <- sse/n_p
Transpuesta <- t(M_X)
A <- solve(Transpuesta%*%M_X)
resultado <- (mse[1]*A)^(1/2)
colnames(resultado) <- c(1,2,3,4)
row.names(resultado) <- c(1,2,3,4)
return (resultado)
}Funcion que calcula el S_Gorro_Yh
S_Gorro_Yh <- function (M_X,M_Y,M_Xh){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
sse <- SSE(M_X,M_Y)
mse <- sse/n_p
Transpuesta <- t(M_X)
A <- solve(Transpuesta%*%M_X)
resultado_a <- mse[1]*A
Transpuesta_M_Xh <- t(M_Xh)
return ((Transpuesta_M_Xh%*%resultado_a%*%M_Xh)^(1/2))
}Funcion que calcula el valor de Bonferroni
#EL VALOR DE PROBAILIDAD INGRESAR ES ALPHA/2, PARAMETROS TIENE EL VALOR
#DEL NUMERO DE PARAMETROS A CALCULAR LIMITE
Valor_Bonferroni <- function (M_X,M_Y,Probabilidad,Parametros){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
t <- Valor_t(Probabilidad*(1/Parametros),n_p)
}Funcion que determina el intervalo de confianza Bonferroni
#INTERVALO DE CONFIANZA CONJUNTO DE BONFERRONI
#EL VALOR DEL ESTIMADO i ES EL B_i CORRESPONDIENTE, PARA B_1 ES 1
#EL VALOR DE PROBABILIDAD ES ALPHA/2
Intervalo_Confianza_Conjunta_Bonferroni <- function (M_X,M_Y,Probabilidad,Parametros,Estimado_i){
B <- Valor_Bonferroni(M_X,M_Y,Probabilidad,Parametros)
b_i <- Matriz_b(M_X,M_Y)[Estimado_i+1]
s_b <- S_b(M_X,M_Y)[Estimado_i+1,Estimado_i+1]
Izquierda <- b_i - B*s_b
Derecha <- b_i + B*s_b
Resultado <- c(Izquierda,"<Betta<",Derecha)
return (Resultado)
}Funcion que determina el intervalo de confianza de E_Y_h
Intervalo_Confianza_E_Yh <- function (M_X,M_Y,M_Xh,Probabilidad){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
Gorro_Yh <- Matriz_Gorro_Yh(M_X,M_Y,M_Xh)[1]
t <- Valor_t(Probabilidad,n_p)
s_gorro_yh <- S_Gorro_Yh(M_X,M_Y,M_Xh)[1]
Izquierda <- Gorro_Yh - t*s_gorro_yh
Derecha <- Gorro_Yh + t*s_gorro_yh
Resultado <- c(Izquierda,"<E{Yh}<",Derecha)
return (Resultado)
}Funcion que realiza la prueba de hipotesis para B_k
Prueba_Hipotesis_Bk <- function (M_X,M_Y,Probabilidad,Estimado_i){
n <- length(M_Y)
p_1 <- ncol(M_X)-1
n_p <- (n-1)-p_1
M1 <- "Rechazamos Ho"
M2 <- "No rechazamos Ho"
Critico <- Valor_t(Probabilidad,n_p)
b_i <- Matriz_b(M_X,M_Y)[Estimado_i+1]
s_b <- S_b(M_X,M_Y)[Estimado_i+1,Estimado_i+1]
Estadistico_Prueba <- b_i/s_b
print("Estadistico Prueba")
print(Valor_Absoluto(Estadistico_Prueba))
print("Critico")
print(Critico)
if (Valor_Absoluto(Estadistico_Prueba)>Critico){
return (M1)
}else{
return (M2)
}
}Funcion que realiza la prueba de hipotesis para la inferencia de una variable en un modelo
Prueba_Hipotesis_Inferencia_Variable<-function (Modelo_Original_X,Modelo_Nuevo_X,M_Y,Probabilidad){
M1 <- "Las variables describen bien el experimento"
M2 <- "Las variables no describen bien el experimento"
sse_modelo_original <- as.numeric(ANOVA(Modelo_Original_X,M_Y,Probabilidad)[2,2])
sse_modelo_nuevo <- as.numeric(ANOVA(Modelo_Nuevo_X,M_Y,Probabilidad)[2,2])
gl_modelo_original <- as.numeric(ANOVA(Modelo_Original_X,M_Y,Probabilidad)[2,3])
gl_modelo_nuevo <- as.numeric(ANOVA(Modelo_Nuevo_X,M_Y,Probabilidad)[2,3])
A <- (sse_modelo_original-sse_modelo_nuevo)/(gl_modelo_original-gl_modelo_nuevo)
B <- (sse_modelo_nuevo/gl_modelo_nuevo)
Estadistico_Prueba <- A/B
Critico <- Valor_F(Probabilidad,gl_modelo_original-gl_modelo_nuevo,gl_modelo_nuevo)
print("Estadistico de Prueba")
print(Estadistico_Prueba)
print("Valor Critico")
print(Critico)
if (Estadistico_Prueba>Critico){
return(M1)
}else{
return(M2)
}
}Añadimos un ejemplo ilustrativo para el uso de las funciones anterior definidas.
Suponga que la base de datos que tenemos es la siguiente.
D_1 <- read.csv("Datos_Problema_1.csv")
D_1## Equipo y x1 x2 x3 x4 x5 x6 x7 x8 x9
## 1 1 10 2113 1985 38.9 64.7 4 868 59.7 2205 1917
## 2 2 11 2003 2855 38.8 61.3 3 615 55.0 2096 1575
## 3 3 11 2957 1737 40.1 60.0 14 914 65.6 1847 2175
## 4 4 13 2285 2905 41.6 45.3 -4 957 61.4 1903 2476
## 5 5 10 2971 1666 39.2 53.8 15 836 66.1 1457 1866
## 6 6 11 2309 2927 39.7 74.1 8 786 61.0 1848 2339
## 7 7 10 2528 2341 38.1 65.4 12 754 66.1 1564 2092
## 8 8 11 2147 2737 37.0 78.3 -1 797 58.9 2476 2254
## 9 9 4 1689 1414 42.1 47.6 -3 714 57.0 2577 2001
## 10 10 2 2566 1838 42.3 54.2 -1 797 58.9 2476 2254
## 11 11 7 2363 1480 37.3 48.0 19 984 68.5 1984 2217
## 12 12 10 2109 2191 39.5 51.9 6 819 59.2 1901 1686
## 13 13 9 2295 2229 37.4 53.6 -5 1037 58.8 1761 2032
## 14 14 9 1932 2204 35.1 71.4 3 986 58.6 1709 2025
## 15 15 6 2213 2140 38.8 58.3 6 819 59.2 1901 2686
## 16 16 5 1722 1730 36.6 52.6 -19 791 54.4 2288 2835
## 17 17 5 1498 2072 35.3 59.3 -5 776 49.6 2072 1914
## 18 18 5 1873 2929 41.1 55.3 10 789 54.3 2861 2496
## 19 19 6 2118 2268 38.6 69.6 6 582 58.7 2411 2670
## 20 20 4 1775 1983 39.3 78.3 7 901 51.7 2289 2202
## 21 21 3 1904 1792 39.7 38.1 -9 734 61.9 2203 1988
## 22 22 3 1929 1606 39.7 68.8 -21 627 52.7 2592 2324
## 23 23 4 2080 1492 35.5 68.8 -8 722 57.8 2053 2550
## 24 24 10 2301 2835 35.3 74.1 2 683 59.7 1979 2110
## 25 25 6 2040 2416 38.7 50.0 0 576 54.9 2048 2628
## 26 26 8 2447 1638 39.9 57.1 -8 848 65.3 1786 1776
## 27 27 2 1416 2649 37.4 56.3 -22 684 43.8 2876 2425
## 28 28 0 1503 1503 39.3 47.0 -9 875 53.5 2560 2241Añadimos un valor de Alpha de 0.05
Alpha <- 0.05Definimos las variables
x1 <- D_1$x1
x2 <- D_1$x2
x3 <- D_1$x3
x4 <- D_1$x4
x5 <- D_1$x5
x6 <- D_1$x6
x7 <- D_1$x7
x8 <- D_1$x8
x9 <- D_1$x9
y <- D_1$y
Matriz_Y <- cbind(y)Recordemos que definimos una función M_X, esta funcion la modificamos para introducir distintos modelos y analizarlos.
#CUATRO VARIABLE
M_X_4 <- function (Vector_X1,Vector_X2,Vector_X3,Vector_X4){
Vector <- c()
for (i in 1:length(Vector_X1)){
Vector <- c(Vector,1)
}
M <- cbind(Vector,Vector_X1,Vector_X2,Vector_X3,Vector_X4)
return (M)
}
#OCHO VARIABLES
M_X_8 <- function (Vector_X1,Vector_X2,Vector_X3,Vector_X4,Vector_X5,Vector_X6,Vector_X7,Vector_X8){
Vector <- c()
for (i in 1:length(Vector_X1)){
Vector <- c(Vector,1)
}
M <- cbind(Vector,Vector_X1,Vector_X2,Vector_X3,Vector_X4,Vector_X5,Vector_X6,Vector_X7,Vector_X8)
return (M)
}
#NUEVE VARIABLES
M_X_9 <- function (Vector_X1,Vector_X2,Vector_X3,Vector_X4,Vector_X5,Vector_X6,Vector_X7,Vector_X8,Vector_X9){
Vector <- c()
for (i in 1:length(Vector_X1)){
Vector <- c(Vector,1)
}
M <- cbind(Vector,Vector_X1,Vector_X2,Vector_X3,Vector_X4,Vector_X5,Vector_X6,Vector_X7,Vector_X8,Vector_X9)
return (M)
}El modelo completo es definido de la siguiente manera.
Modelo_Nuevo_9 <- M_X_9(x1,x2,x3,x4,x5,x6,x7,x8,x9)Para determinar la recta de regresión ajustada a los datos calculamos los valores de los b_i
Matriz_b9 <- function(M_X,M_Y){
T_M_X <- t(M_X)
A1 <- T_M_X%*%M_X
A <- solve(A1)
B <- T_M_X%*%M_Y
Resultado <- A%*%B
row.names(Resultado) <- c("b0","b1","b2","b3","b4","b5","b6","b7","b8","b9")
colnames(Resultado) <- c("Coeficientes de Regresion")
return (Resultado)
}
Matriz_b9(Modelo_Nuevo_9,Matriz_Y)## Coeficientes de Regresion
## b0 -1.283255e+01
## b1 8.168413e-04
## b2 3.962530e-03
## b3 1.752412e-02
## b4 4.813196e-02
## b5 -2.738455e-02
## b6 3.038332e-03
## b7 2.389176e-01
## b8 -2.811323e-03
## b9 -1.922779e-03Calculamos el ANOVA del modelo completo
ANOVA(Modelo_Nuevo_9,Matriz_Y,Alpha)## FUENTE_VARIACION SS GL MS
## [1,] "Regresion" "265.961898807713" "9" "29.5513220897459"
## [2,] "Error" "61.0023869065727" "18" "3.3890214948096"
## [3,] "Total" "326.964285714286" "27" "-"
## ESTADISTICO_PRUEBA VALOR_CRITICO
## [1,] "8.71972105665449" "2.45628114915927"
## [2,] "-" "-"
## [3,] "-" "-"Analizamos la contribucion individual de cada variable al modelo completo.
#PARA x1
Contribucion(M_X_8(x2,x3,x4,x5,x6,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 0.5719314#PARA x2
Contribucion(M_X_8(x1,x3,x4,x5,x6,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 81.83434#PARA x3
Contribucion(M_X_8(x1,x2,x4,x5,x6,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 0.01606903#PARA x4
Contribucion(M_X_8(x1,x2,x3,x5,x6,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 4.261248#PARA x5
Contribucion(M_X_8(x1,x2,x3,x4,x6,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 1.091963#PARA x6
Contribucion(M_X_8(x1,x2,x3,x4,x5,x7,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 2.853179#PARA x7
Contribucion(M_X_8(x1,x2,x3,x4,x5,x6,x8,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 9.210075#PARA x8
Contribucion(M_X_8(x1,x2,x3,x4,x5,x6,x7,x9),Modelo_Nuevo_9,Matriz_Y)## [1] 8.751425#PARA x9
Contribucion(M_X_8(x1,x2,x3,x4,x5,x6,x7,x8),Modelo_Nuevo_9,Matriz_Y)## [1] 8.198293