10 Nomograma

Esta é uma pré-visualização de arquivo. Entre para ver o arquivo original

########################################
#
# CONTROLE ESTATÍSTICO DA QUALIDADE
# NOMOGRAMA
#
# Prof. Dr. Vinícius Quintas
# Departamento de Estatística - UFPE
#
########################################
NQA = 0.02
P_NQA = 0.80
NQR = 0.06
P_NQR = 0.20
#--------------------
# Amostragem Única 
#--------------------
lim_c = 50
lim_n = 1000
Pa_NQA = matrix(0,nrow=lim_n, ncol=lim_c+1)
Pa_NQR = matrix(0,nrow=lim_n, ncol=lim_c+1)
for(n in 1:lim_n)
{ for(c in 0:lim_c)
 { if(c<=n)
 { Pa_NQA[n,c+1] = pbinom(c,n,NQA)
 Pa_NQR[n,c+1] = pbinom(c,n,NQR)
 }
 } 
}
epsilon=seq(0.001,0.05,0.001)
for(i in 1:length(epsilon))
{ dat_NQA = which(abs(Pa_NQA - P_NQA) < epsilon[i])
 dat_NQR = which(abs(Pa_NQR - P_NQR) < epsilon[i])
 inter = intersect(dat_NQA,dat_NQR)
 if(sum(inter)!=0)
 { values_c = ceiling(inter/lim_n)-1
 values_n = inter - floor(inter/lim_n)*lim_n
 print(matrix(nrow=3, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Planejamento Amostral")),
 paste(c("NQA ",sprintf("%0.2f", NQA)," - P(NQA) ",sprintf("%0.2f",P_NQA)),collapse=""),
 paste(c("NQR ",sprintf("%0.2f", NQR)," - P(NQR) ",sprintf("%0.2f",P_NQR)),collapse="")))
 ))
 print(matrix(nrow=6, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Plano Amostral")),
 paste(c("n = ",sprintf("%0.0f", values_n)),collapse=""),
 paste(c("c = ",sprintf("%0.0f", values_c)),collapse=""),
 paste(c("PA(NQA) = ",sprintf("%0.3f", pbinom(values_c,values_n,NQA)),
 " ; PA(NQR) = ",sprintf("%0.3f",pbinom(values_c,values_n,NQR))),
 collapse=""),
 paste(c("|P(NQA)-PA(NQA)| < ",sprintf("%0.3f", epsilon[i])),collapse=""),
 paste(c("|P(NQR)-PA(NQR)| < ",sprintf("%0.3f", epsilon[i])),collapse="")))))
 { 
 epsilon = epsilon[i]*1.5
 dat_NQA = which(abs(Pa_NQA - P_NQA) < epsilon)
 dat_NQR = which(abs(Pa_NQR - P_NQR) < epsilon)
 inter = intersect(dat_NQA,dat_NQR)
 values_c = ceiling(inter/lim_n)-1
 values_n = inter - floor(inter/lim_n)*lim_n
 print(matrix(nrow=4, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Outros Planos Amostrais")),
 paste(c("|P(NQA)-PA(NQA)| < ",sprintf("%0.3f", epsilon)),collapse=""),
 paste(c("|P(NQR)-PA(NQR)| < ",sprintf("%0.3f", epsilon)),collapse=""),""))))
 print(data.frame(list("n" = values_n,
 "c" = values_c,
 "PA_NQA" = round(pbinom(values_c,values_n,NQA),digits = 3), 
 "P_NQA" = P_NQA,
 "PA_NQR" = round(pbinom(values_c,values_n,NQR),digits = 3),
 "P_NQR" = P_NQR )))
 }
 break
 }
}
 
#--------------------
# Amostragem Dupla 
#--------------------
lim_c1 = 4
lim_c2 = 4
lim_n1 = 70
lim_n2 = 110
Pa1_NQA = matrix(0,nrow=lim_n1, ncol=lim_c1+1)
Pa1_NQR = matrix(0,nrow=lim_n1, ncol=lim_c1+1)
Pa2_NQA = array(0,dim=c(lim_n1,lim_c1+1,lim_n2,lim_c2+1))
Pa2_NQR = array(0,dim=c(lim_n1,lim_c1+1,lim_n2,lim_c2+1))
Pa_NQA = array(0,dim=c(lim_n1,lim_c1+1,lim_n2,lim_c2+1))
Pa_NQR = array(0,dim=c(lim_n1,lim_c1+1,lim_n2,lim_c2+1))
for(n1 in 1:lim_n1)
{ for(c1 in 0:lim_c1)
 { if(c1<=n1)
 { Pa1_NQA[n1,c1+1] = pbinom(c1,n1,NQA)
 Pa1_NQR[n1,c1+1] = pbinom(c1,n1,NQR)
 }
 }
}
for(n1 in 1:lim_n1)
{ for(c1 in 0:lim_c1)
 { for(n2 in 1:lim_n2)
 { for(c2 in 0:lim_c2)
 { if((c1<=n1)&&(c2<=(n1+n2))&&(c2>c1))
 { sum1 = sum2 = 0
 for(d1 in (c1+1):c2)
 { for(d2 in 0:(c2-d1))
 { aux_1 = dbinom(d1,n1,NQA)*dbinom(d2,n2,NQA)
 sum1 = sum1 + aux_1
 aux_2 = dbinom(d1,n1,NQR)*dbinom(d2,n2,NQR)
 sum2 = sum2 + aux_2 
 }
 }
 Pa2_NQA[n1,c1+1,n2,c2+1] = sum1
 Pa2_NQR[n1,c1+1,n2,c2+1] = sum2
 } 
 }
 }
 } 
} 
for(n2 in 1:lim_n2)
{ for(c2 in 0:lim_c2)
 { if(sum(Pa2_NQA[,,n2,c2+1])!=0)
 { Pa_NQA[,,n2,c2+1] = Pa2_NQA[,,n2,c2+1] + Pa1_NQA
 }
 if(sum(Pa2_NQR[,,n2,c2+1])!=0)
 { Pa_NQR[,,n2,c2+1] = Pa2_NQR[,,n2,c2+1] + Pa1_NQR
 }
 }
}
#Pa_NQA
epsilon=seq(0.001,0.05,0.001)
for(i in 1:length(epsilon))
{ dat_NQA = which(abs(Pa_NQA - P_NQA) < epsilon[i])
 dat_NQR = which(abs(Pa_NQR - P_NQR) < epsilon[i])
 inter = intersect(dat_NQA,dat_NQR)
 if(sum(inter)!=0)
 { bloco = ceiling(inter/(lim_n1*(lim_c1+1)))
 linha = inter - (bloco-1)*lim_n1*(lim_c1+1)
 n1 = linha - lim_n1*(ceiling(linha/lim_n1)-1)
 c1 = ceiling(linha/lim_n1)-1
 n2 = bloco - (ceiling(bloco/lim_n2)-1)*lim_n2
 c2 = floor(bloco/lim_n2)
 print(matrix(nrow=3, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Planejamento Amostral")),
 paste(c("NQA ",sprintf("%0.2f", NQA)," - P(NQA) ",sprintf("%0.2f",P_NQA)),collapse=""),
 paste(c("NQR ",sprintf("%0.2f", NQR)," - P(NQR) ",sprintf("%0.2f",P_NQR)),collapse="")))
 ))
 print(matrix(nrow=8, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Plano Amostral")),
 paste(c("n1 = ",sprintf("%0.0f", n1)),collapse=""),
 paste(c("c1 = ",sprintf("%0.0f", c1)),collapse=""),
 paste(c("n2 = ",sprintf("%0.0f", n2)),collapse=""),
 paste(c("c2 = ",sprintf("%0.0f", c2)),collapse=""),
 paste(c("PA(NQA) = ",sprintf("%0.3f", Pa_NQA[n1,c1+1,n2,c2+1]),
 " ; PA(NQR) = ",sprintf("%0.3f",Pa_NQR[n1,c1+1,n2,c2+1])),
 collapse=""),
 paste(c("|P(NQA)-PA(NQA)| < ",sprintf("%0.3f", epsilon[i])),collapse=""),
 paste(c("|P(NQR)-PA(NQR)| < ",sprintf("%0.3f", epsilon[i])),collapse="")))))
 { 
 epsilon = epsilon[i]*1.5
 dat_NQA = which(abs(Pa_NQA - P_NQA) < epsilon)
 dat_NQR = which(abs(Pa_NQR - P_NQR) < epsilon)
 inter = intersect(dat_NQA,dat_NQR)
 bloco = ceiling(inter/(lim_n1*(lim_c1+1)))
 linha = inter - (bloco-1)*lim_n1*(lim_c1+1)
 n1 = linha - lim_n1*(ceiling(linha/lim_n1)-1)
 c1 = ceiling(linha/lim_n1)-1
 n2 = bloco - (ceiling(bloco/lim_n2)-1)*lim_n2
 c2 = floor(bloco/lim_n2)
 
 aux1=aux2=NULL
 for(j in 1:length(c2))
 { aux1[j] = Pa_NQA[n1[j],(c1[j]+1),n2[j],c2[j]+1]
 aux2[j] = Pa_NQR[n1[j],(c1[j]+1),n2[j],c2[j]+1]
 } 
 print(matrix(nrow=4, ncol=0, byrow=TRUE,dimnames =list(c(
 paste(c("Outros Planos Amostrais")),
 paste(c("|P(NQA)-PA(NQA)| < ",sprintf("%0.3f", epsilon)),collapse=""),
 paste(c("|P(NQR)-PA(NQR)| < ",sprintf("%0.3f", epsilon)),collapse=""),""))))
 print(data.frame(list("n1" = n1,"c1" = c1, "n2" = n2, "c2" = c2,
 "PA_NQA" = round(aux1,digits=3), "P_NQA" = P_NQA, "PA_NQR" = round(aux2,digits=3), "P_NQR" = P_NQR )))
 }
 break
 }
}
#EXEMPLO COMPARATIVO
N=73
C=2
n1=50
n2=52
c1=1
c2=2
P2A=0
p=seq(0,0.12,length=130)
p=0.02
p=0.06
P1A = pbinom(c1,n1,p)
P1R = 1 - pbinom(c2,n1,p)
for(d1 in (c1+1):c2)
{
 for(d2 in 0:(c2-d1))
 {
 aux = dbinom(d1,n1,p)*dbinom(d2,n2,p)
 P2A = P2A + aux
 }
} 
PA = P1A + P2A
plot(p,P1A,type="l",ylim=c(0,1),ylab="Prob(Aceitação)",xlab="p")
par(new=TRUE)
plot(p,pbinom(2,73,p),type="l",ylim=c(0,1),lty=2,lwd=2,ylab="",xlab="")
grid()