# ***********************************************************************************
# MQCI: "MULTIVARIATE QUALITY CONTROL IMPLEMENTATION"
# Macro: mqci.ssc
# Versão: 1.0 
# Ano: 2005
# Autores: Fernanda Karine Ruiz Colenghi (Bacharelanda em Estatística)
#          Sueli Aparecida Mingoti (Profa. Adjunta- UFMG - orientadora)

# Descrição do programa é a seguinte: Este programa implementa os testes T2 
# de Hotelling e o de Hayter e Tsui para controle de processos multivariados, 
# tanto para observações individuais como para subgrupos racionais.
# Permite também que o usuário identifique as variáveis causadoras da
# falta de controle do processo. 
# Este programa faz análise de dados para as fases I e II de controle de qualidade.

# Nota: Este programa funciona em S-Plus versões 2000 ou superior.
# ***********************************************************************************


######################## Funcao que calcula as estatisticas basicas ##################################

FEstatisticasBasicas_function(dados1,m,n,p,valores,dadmedia,dadcov)
{
	dados2_dados1[ ,-1] 
	dadm_data.matrix(dadmedia)
	dadc_data.matrix(dadcov)
 
	if (n>1) 
	{
		if (valores==1) #usar a estimacao
		{
			# calculo de xjkbarra
			x1_aggregate.data.frame(dados1,dados1[ ,1],mean)
			x2_x1[ ,-c(1, 2)] 
			xbarra_data.matrix(x2) # matriz m x p
			xbarra_t(xbarra) # matriz p x m
		
			# calculo de xbarrabarra
			med_apply(dados2,2,mean) # nao ve como matriz
			medi_t(med) # matriz 1 x p (comando interno)
			media_t(medi) # matriz p x 1

			# calculo de S e R 
			dadosn_list
			dadosr_list
			for (k in 1:m)
			{
				vardados_dados1[dados1[ ,1]==k, ]
				vardados_vardados[ ,-c(1)]
				covariancia_var(vardados)
				dadosn[[k]]_covariancia
				correlacao_cor(vardados)
				dadosr[[k]]_correlacao
			}
			somas_0
			somar_0
			for(k in 1:m)
				{
					somas_somas+dadosn[[k]]
					somar_somar+dadosr[[k]]
				}
			s_somas/m  # matriz p x p
			corr_somar/m # matriz p x p
			# matriz precisa ser inversivel para funcionar????
			
			# vetor com variancias
			diagonal_diag(s) # nao ve como matriz
			varian_t(diagonal) # matriz 1 x p (comando interno)
			varianc_t(varian) # matriz p x 1

		}
		else # (valores==2)
		{                  
			# calculo de xjkbarra
			x1_aggregate.data.frame(dados1,dados1[ ,1],mean)
			x2_x1[ ,-c(1,2)] 
			xbarra_data.matrix(x2) # matriz m x p
			xbarra_t(xbarra) # matriz p x m
		
			# calculo de xbarrabarra
			media_data.matrix(dadmedia) # matriz 1 x p
			media_t(media) # matriz p x 1 

			# calculo de S
			s_dadc # matriz p x p
								
			# vetor com variancias
			diagonal_diag(s) # nao ve como matriz
			VetDiagInversa_1/diagonal
			DiagID_diag(VetDiagInversa)# transforma numa matriz diagonal
			
			# calculo de R
			corr_DiagID%*%s%*%DiagID
			varian_t(diagonal) # matriz 1 x p (comando interno)
			varianc_t(varian) # matriz p x 1

		}
	}
	else #n=1		
	{
		if (valores==1)
		{
			# calculo de xjkbarra que e o valor observado
			xbarra_data.matrix(dados2) # matriz m x p
			xbarra_t(xbarra) # matriz p x m		
		
			# calculo de xbarrabarra
			med_apply(dados2,2,mean) # nao ve como matriz
			medi_t(med) # matriz 1 x p (comando interno)
			media_t(medi) # matriz p x 1
					
			# vetor com variancias
			varia_apply(dados2,2,var) # nao ve como matriz
			varian_t(varia) # matriz 1 x p (comando interno)
			varianc_t(varian) # matriz p x 1

			# calculo da matriz de covariancia amostral e correlacao
			s_var(dados2) # matriz p x p
			corr_cor(dados2)	
		}
		else # valores==2
		{
			# calculo de xjkbarra que e o valor observado
			xbarra_data.matrix(dados2) # matriz m x p
			xbarra_t(xbarra) # matriz p x m

			# calculo de xbarrabarra
			media_dadm # matriz 1 x p
			media_t(media) # matriz p x 1 
				
			# calculo de S e R
			s_dadc # matriz p x p
								
			# vetor com variancias
			diagonal_diag(s) # nao ve como matriz
			VetDiagInversa_1/diagonal
			DiagID_diag(VetDiagInversa)
			
			# calculo de R
			corr_DiagID%*%s%*%DiagID
			varian_t(diagonal) # matriz 1 x p (comando interno)
			varianc_t(varian) # matriz p x 1
		}
	}
	return(media,xbarra,s,corr,varianc)
}

############################# Funcao que calcula os Tk's para cada amostra #############################

FCalculaTk_function(dados1,m,n,p,valores,dadmedia,dadcov)
{
	resutados_FEstatisticasBasicas(dados1,m,n,p,valores,dadmedia,dadcov)
	media_resutados$media
	s_resutados$s
	xbarra_resutados$xbarra
	dif_matrix(0, p, m)
	T2_rep(0, m)
	for(k in 1:m) 
	{
		dif[,k]_xbarra[, k]-media # matriz p x 1
		diftrans_ # matriz 1 x p
		inversa_solve(s) # martriz p x p
		T2[k]_n*(t(dif[, k])%*%inversa%*% dif[, k]) 
		# matriz(1xp)*(pxp)*(px1)=matriz(1x1)=numero
	}
	return(T2)
}

############################# Funcao que calcula Cralfa ################################################
Fcralfa_function(dados1,alfa,p,corr)
{
	v_rmvnorm(10000,mean=rep(0,p),cov=corr)
	g_abs(v)
	h_apply(g,1,max)
	o_sort(h)
	conf_10000*(1-alfa)
	cralfa_o[conf]
	return(cralfa)
}

#####################################################################################################
 

################################### Funcao que faz a padronizacao ####################################

FPadroniza_function(dados1,m,n,p,media,corr,varianc)
{
	col_dados1[ ,1]
	dados2_dados1[ ,-1] 
	desvi_sqrt(varianc) # p*1
	
	tmedia_t(media) # 1*p
	tot_n*m 
	tmedi_rep(tmedia,tot)
	tmedia_matrix(tmedi,ncol=p,byrow=T) # p*(mxn)
	tdesvi_t(desvi) # 1*p
	tdesv_rep(tdesvi,tot)
	tdesvi_matrix(tdesv,ncol=p,byrow=T) # p*(mxn)
	
	pad_dados2-tmedia # p*(mxn)
	pad_abs(pad) # modulo dos valores
	padr_pad/tdesvi
	dadpadr_cbind(col,padr) # (p+1)*(mxn)
	return(dadpadr) # dados padronizados
}
	
######################################################################################################


################################ inicio do programa principal #########################################

controlem_function(dados,i,m,n,p,alfa,fase,valores,dadmedia,dadcov,estcralfa,cr,saida)
{
# dados do data.sheet
# i coluna de identificacao das amostras (subgrupos)
# m e o numero de amostras
# n e o tamanho das amostras(subgrupos)
# p e o numero de variaveis
# alfa e o nivel de significancia 
# fase= 1 para fazer o teste na fase 1 e 2 se o usuario quer o teste apenas na fase 2
# Valores= 1 estimar: dadmedia=0, dadcov=0 
# Valores= 2 buscar o banco de dados em dadmedia e dadcov
# o programa funciona com os dados representados da forma tal que ha 
# a identificacao da amostra na primeira coluna, nas outras estao as p variaveis
# estcralfa= 1 estimar com os dados: cr=0
# estcralfa= 2 usar o valor de cralfa fornecido pelo usuário em cr
# e relevante mencionar que o valor de cralfa so pode ser fornecido pelo usuário na fase 2
# se o usuário deseja o resumo das saidas ele coloca 1 em saida
# se o usuário deseja todas saídas ele coloca 2 em saida
	dados1_data.matrix(dados)
	dadm_data.matrix(dadmedia)
	dadc_data.matrix(dadcov)

	dados2_dados1[ ,-c(i)]
	id_dados1[,i]
	dados1_cbind(id, dados2)
	problema2_0
	eme_m
	
# Autoria do programa
	cat("\n")
	cat("*****************************************************************************************", fill=T)
	cat(" MQCI: 'MULTIVARIATE QUALITY CONTROL IMPLEMENTATION' ", fill=T)
	cat(" Macro: mqci.ssc", fill=T)
	cat(" Versão: 1.0", fill=T)
	cat(" Ano: 2005", fill=T) 
	cat(" Autores: Fernanda Karine Ruiz Colenghi (Bacharelanda em Estatistica) ", fill=T) 
	cat("          Sueli Aparecida Mingoti (Profa. Adjunta- UFMG - orientadora)", fill=T) 
	cat(" Descricao do programa e a seguinte: Este programa implementa os testes T2 de Hotelling ", fill=T) 
	cat(" e o de Hayter e Tsui para controle de processos multivariados, tanto para observacoes " ,fill=T)   
	cat(" individuais como para subgrupos racionais. Permite tambem que o usuario identifique ", fill=T)  
	cat(" 	as variaveis causadoras da falta de controle do processo.", fill=T) 
	cat(" Este programa faz analise de dados para as fases I e II de controle de qualidade.", fill=T) 
	cat(" Nota: Este programa funciona em S-Plus versões 2000 ou superior.", fill=T)
	cat("*****************************************************************************************", fill=T)
	cat("\n")

############ Procedimento 1: calcula as estatisticas basicas de media e de variancia ################

	resutados_FEstatisticasBasicas(dados1,m,n,p,valores,dadmedia,dadcov)
	media_resutados$media
	s_resutados$s
	xbarra_resutados$xbarra
	corr_resutados$corr 
	varianc_resutados$varianc

	
	# imprime o vetor de medias e a matriz de covariancias
	if (valores==1)
	{
		cat("A media amostral observada e dada por: ",fill=T)
		print(xbarra)
		cat("O vetor de medias estimado e dado por: ",fill=T)
		print(media)
		cat("A matriz de covariancias estimada e dada por: ",fill=T)
		print(s)
		cat("A matriz de correlacao estimada e dada por: ",fill=T)
		print(corr)

	}
	else # valores=2
	{
		cat("A media amostral observada e dada por: ",fill=T)
		print(xbarra)
		cat("O vetor de medias fornecido e: ",fill=T)
		print(media)
		cat("A matriz de covariancias fornecida e: ",fill=T)
		print(s)
		cat("A matriz de correlacao calculdada e dada por: ",fill=T)
		print(corr)
	}

######################################### Fim do procedimento 1 ########################################

	cat("******************************************************************************************", fill=T)
	cat("\n")

############### Procedimento 2: Calculo das estatísticas Tks para cada amostra ###########################

	T2_FCalculaTk(dados1,m,n,p,valores,dadmedia,dadcov)
	
	if(fase==1&valores==1|fase==2)
	{
		if (saida==2)
		{
			cat("Valores de T2: ",fill=T)
			print(T2)
		}
	}
############################### Fim do procedimento 2 ###################################################
	if (saida==2)
	{
		cat("****************************************************************************************", fill=T)
		cat("\n")
	}

################# Procedimento 3: Calculo do Cralfa ########################################################
	if (fase==1)
	{
		cralfa_Fcralfa(dados1,alfa,p,corr)
		if(estcralfa==1)
		{
			cat("Valor de cralfa: ",cralfa,fill=T)
		}
		else # estcralfa==2
		{
			cat(" Cralfa nao pode ser fornecido. O programa usara a estimativa de cralfa: ",cralfa)
			cat(" obtido com os dados", fill=T)
		}
	}
	else # fase==2
	{
		if(estcralfa==1)
		{
			cralfa_Fcralfa(dados1,alfa,p,corr)
			cat("Valor de cralfa: ")
			print(cralfa)
			cat("\n")
		}
		else # estcralfa==2
		{
			cralfa_cr
			cat("Valor de cralfa: ",cralfa," fornecido pelo usuario",fill=T)
		}
	}

############################### Fim do procedimento 3 ###################################################

	cat("******************************************************************************************", fill=T)
	cat("\n")

###################### Procedimento 4: gráficos de controle da fase 1 ###############################
	if (fase==1)
	{
		if (valores==1)
		{
			if(n>1)
			{
				f_qf(1-alfa,p,m*n-m-p-1)
				LSC1_p*(m-1)*(n-1)*f/(m*n-m-p+1)
			}
			else # n=1
			{
				b_qbeta(1-alfa,p/2,(m-p-1)/2) # conferir
				LSC1_(m-1)*(m-1)*b/m
			}
			cat("Limite superior de controle de T2 de Hotelling e: ",LSC1,fill=T)
			retiradas_0
			problema1_0
			probf1_rep(0, m)
			if (max(T2)<=LSC1)
			{
				# não existe nenhum T2 fora dos limites de controle
				cat("O processo esta sob controle estatistico sob o teste T2 de Hotelling.\n")
			}
			else
			{
				cat("Existe(m) amostra(s) fora dos limites de controle por T2 de Hotelling", fill=T)
				cat(" As amostras fora dos limites de controle por Hotelling sao: \n")
			}
			while(max(T2)>LSC1)
			{
				k_m
				chama_1
				while(k>=1)
				{
					if(	T2[k]>LSC1)
					{
						ExcluiS_k*n
						ExcluiI_ExcluiS-n+1
						dados1_dados1[-c(ExcluiI:ExcluiS),]
						T2_T2[-k]
						m_m-1
						probf1[k]_1
						print(k)
						retiradas_retiradas+1
					}
					k_k-1
				}
				#cat("Dados sem as amostras problematicas: ", fill=T)
				#print(dados1)
				k_m
				while(k>=1)
				{
					k1_0
					while(k1<n)
					{
						Pos_(k*n)-k1
						dados1[Pos,1]_k
						k1_k1+1
					}
					k_k-1
				}
				problema1_max(probf1)
				#cat("\n dados1: \n")
				# calculo das novas estatisticas basicas
				# fazer novo calculo das estatisticas basicas 
				#cat("Novo banco de dados para analise: ",fill=T)
				#print(dados1)

				resutados_FEstatisticasBasicas(dados1,m,n,p,valores,dadmedia,dadcov)
 				media_resutados$media
				s_resutados$s
				xbarra_resutados$xbarra
				corr_resutados$corr
				
				# Calculo dos novos T2
				#fazer novo calculo dos T2
				T2_FCalculaTk(dados1,m,n,p,valores,dadmedia,dadcov)
				#cat("\n T2: ",T2,fill=T)
				
				# escrever a nova funcao
				if (n>1)
				{
					f_qf(1-alfa,p,m*n-m-p-1)
					LSC1_p*(m-1)*(n-1)*f/(m*n-m-p+1)
				}
				else # n=1
				{
					b_qbeta(1-alfa,p/2,(m-p-1)/2) # conferir
					LSC1_(m-1)*(m-1)/m
				}
			}
			
			# Grafico de controle
			LSC1Vet_rep(LSC1,m)  
			tsplot(T2,LSC1Vet,xlab="amostra",ylab="T2")
			title("Grafico de Controle T2 de Hotelling para a Fase 1")
			
			if(chama==1)
			{
				cat("\n O numero de amostras em que iremos trabalhar é: ",m, fill=T)
				cat("Limite superior de controle final de T2 de Hotelling : ",LSC1,fill=T)
				cat(" Estatisticas descritivas dos dados apos a implementacao do T2 de Hotelling na fase 1.\n")
				cat("\n media: ",fill=T)
				print(media)
				cat("\n covariancia: ",fill=T)
				print(s)
				cat("\n correlacao: ", fill=T)
				print(corr)
			}

			cat("*************************************************************************************", fill=T)
			cat("\n")
			# Chamar de novo a funcao se for necessario
			dados1_data.matrix(dados)
			dados2_dados1[ ,-c(i)]
			id_dados1[,i]
			dados1_cbind(id, dados2)

			m_eme
			
			#### Algoritmo cralfa #########
			dadpadr_FPadroniza(dados1,m,n,p,media,corr,varianc)
			
			maxi_aggregate.data.frame(dadpadr,dadpadr[ ,1],mean)
			maxi_maxi[ ,-c(1, 2)]  # matriz m x p
			maxim_apply(maxi,1,max) # nao ve como matriz
			maxim_t(maxim) # matriz 1 x m (comando interno)
			maximopd_t(maxim) # matriz m x 1
			
			if (saida==2)
			{
				cat(" O maximo dos dados padronizados em cralfa e: \n")
				print(maximopd)
			}

			if(max(maximopd)<=cralfa)
			{
				# nao existe nenhum valor maior que o cralfa 
				cat("O processo esta sob controle pela metodologia de Hayter e Tsui.\n")
				chame_0
			}
			else
			{
				cat("Existe(m) amostra(s) fora dos limites de controle pela metodologia de Hayter e Tsui.\n")
				cat(" As amostras fora dos limites de controle em Hayter e Tsui sao: \n")
				# existe pelo menos uma observacao acima dos limites
			}
			retiradas_0
			problema1_0
			probf1_rep(0, m)
		
			while (max(maximopd)>cralfa)
			{
				k_m
				chame_1
				while(k>=1)
				{
					if(	maximopd[k]>cralfa)
					{
						ExcluiS_k*n
						ExcluiI_ExcluiS-n+1
						dados1_dados1[-c(ExcluiI:ExcluiS),]
						maximopd_maximopd[-k]
						m_m-1
						probf1[k]_1
						print(k)
						retiradas_retiradas+1
					}
					k_k-1
				}
				#cat("Dados sem as amostras problematicas: ", fill=T)
				#print(dados1)
				k_m
				while(k>=1)
				{
					k1_0
					while(k1<n)
					{
						Pos_(k*n)-k1
						dados1[Pos,1]_k
						k1_k1+1
					}
					k_k-1
				}
				problema1_max(probf1)
				#cat("\n dados1: \n")
				# calculo das novas estatisticas basicas
				# fazer novo calculo das estatisticas basicas 
				#cat("Novo banco de dados para analise: ",fill=T)
				#print(dados1)
				#cat("\n O numero de amostras que iremos trabalhar e: ",m)

				resutados_FEstatisticasBasicas(dados1,m,n,p,valores,dadmedia,dadcov)
 				media_resutados$media
				s_resutados$s
				xbarra_resutados$xbarra
				corr_resutados$corr
				
				# Calculo do Cralfa
				cralfa_Fcralfa(dados1,alfa,p,corr)
				#cat("Novo cralfa: ",cralfa,fill=T)
				# Padronizar os dados
				dadpadr_FPadroniza(dados1,m,n,p,media,corr,varianc)		
				
				maxi_aggregate.data.frame(dadpadr,dadpadr[ ,1],mean)
				maxi_maxi[ ,-c(1, 2)]  # matriz m x p
				maxim_apply(maxi,1,max) # nao ve como matriz
				maxim_t(maxim) # matriz 1 x m (comando interno)
				maximopd_t(maxim) # matriz m x 1
			}
			
			# Grafico de controle do Cralfa
			LSCralfa_rep(cralfa,m)
			tsplot(maximopd,LSCralfa,xlab="amostra",ylab="estatistica M")
			title("Grafico de Controle de Hayter e Tsui para a Fase 1")

			if(chame==1)
			{
				cat(" Novo cralfa: ", cralfa, fill=T)
				cat(" Estatisticas descritivas dos dados apos a metodologia de Hayter e Tsui na fase 1.\n")
				cat("\n media: ",fill=T)
				print(media)
				cat("\n covariancia: ",fill=T)
				print(s)
				cat("\n correlacao: ", fill=T)
				print(corr)
			}
			else
			{
				cat("\n")
			}

			cat("************************************************************************************", fill=T)
			cat("\n")	
		}
		else # valores==2
		{
			cat("Erro",fill=T)
		}
	}
	
	
################################ Fim do procedimento 4 ######################################################
	
	
################################# Procedimento 5: inicio da fase 2 #########################################
	
	else #fase==2
	{
		if (valores==1)
		{
			if (n>1)
			{
				f1_qf(1-alfa,p,m*n-m-p-1)
				LSC_p*(m+1)*(n-1)*f1/(m*n-m-p+1)
			}
			else # n=1
			{
				f1_qf(1-alfa,p,m-p)
				LSC_p*(m+1)*(m-1)*f1/(m*m-m*p)
			}
		}
		else # valores==2
		{
			q2_qchisq(1-alfa,p) 
			LSC_q2
		}
		cat("Limite superior de controle de T2 de Hotelling: ",LSC,fill=T)
		
		# Construcao dos graficos de controle para a fase 2
		LSCVet_rep(LSC,m)
		tsplot(T2,LSCVet,xlab="amostra",ylab="T2")
		title("Grafico de Controle T2 de Hotelling para a Fase 2")
		
		probf2_rep(0, m)
		for(k in 1:m)	
		{
			if(T2[k]>LSC)
			{
				probf2[k]_1
			}
			else # T2[k]<=LSC
			{
				probf2[k]_0
			}
		}
		problema2_max(probf2)
		r_0
		if(problema2==0)
		{
			cat("O processo esta sob controle estatistico pela metodologia de T2 de Hotelling.",fill=T)
		}
		else # existe algum T2 fora dos limites de controle
		{
			cat(" As amostras fora dos limites de controle em T2 de Hotelling sao: ")
			for(k in 1:m)	
			{
				if(	T2[k]>LSC)
				{
					cat(" ",k)
					if(r==0)
					{
						r_r+1
						xbarran_xbarra[,k]
						xbarran_t(xbarran)
						xbarran_t(xbarran)
						T2probn_T2[k]
					}
					else
					{
						r_r+1
						xbarran_cbind(xbarran,xbarra[ ,k]) # matriz(p*r)
						T2probn_cbind(T2probn,T2[k])
					}
				}
			}	
			#cat("\n xbarran \n")
			#print(xbarran)
			cat("\n")
			T2probn_rep(T2probn,p)
			T2probn_matrix(T2probn,ncol=r,byrow=T)# matriz(P*r)
			#cat("T2probn \n")
			#print(T2probn)
			cat("\n")
			cat("O numero de amostras problematicas e: ",r,fill=T)
		}
		cat("****************************************************************************************", fill=T)
		cat("\n")
		
		#### Algoritmo cralfa #########
		dadpadr_FPadroniza(dados1,m,n,p,media,corr,varianc)
		
		maxi_aggregate.data.frame(dadpadr,dadpadr[ ,1],mean)
		maxi_maxi[ ,-c(1, 2)]  # matriz m x p
		maxim_apply(maxi,1,max) # nao ve como matriz
		maxim_t(maxim) # matriz 1 x m (comando interno)
		maximopd_t(maxim) # matriz m x 1
		#cat("Cralfa: ", cralfa, fill=T)
		
		
		# Grafico de controle do Cralfa
		LSCralfa_rep(cralfa,m)
		tsplot(maximopd,LSCralfa,xlab="amostra",ylab="estatistica M")
		title("Grafico de Controle de Hayter e Tsui para a Fase 2")
		if (saida==2)
		{
			cat("Maximo dos dados padronizados em Hayter e Tsui: ",fill=T)
			print(maxim)
			cat("\n")
		}
		ExisteMaior_0
		for(i in 1:m)
		{
			if(maximopd[i]>=cralfa)
				ExisteMaior_1
		}
		
		if(ExisteMaior==0)  ######## olhar maximo do vetor
		{
			# nao existe nenhum valor maior que o cralfa 
			cat("O processo esta sob controle pela metodologia de Hayter e Tsui.\n")
			cat("*************************************************************************************", fill=T)
			cat("\n")

		}
		else # existe pelo menos uma observacao acima dos limites 
		{
			cat(" Existe pelo menos uma obsevacao de forma padronizada que apresenta limites ")
			cat("maiores que o cralfa.",fill=T)
			varianci_t(varianc) # 1*p
			desvio_sqrt(varianci)# 1*p
			medcr_t(media) # 1*p
			LICr_medcr-(cralfa*desvio)
			LSCr_medcr+(cralfa*desvio)
			cat(" Matriz com os limites inferior e superior de controle de Hayter e Tsui ")
			cat("para todas variaveis: ",fill=T)
			tLICr_t(LICr) # p*1
			tLSCr_t(LSCr) # p*1
			limit_cbind(tLICr,tLSCr)# p*2
			print(limit)
			
			# serve tanto para n=1 quanto n>1
			posic_rep(0, m*n)
			
			cat("\n As observacoes fora dos limites de controle em Hayter e Tsui sao: ",fill=T)
			cat(" Variavel,amostra, posicao da observacao na amostra e valor observado", fill=T)
			for(ob in 1:(m*n))	
			{
				for(j in 1:p)
				{
					if(dados1[ob,j+1]<LICr[j]|dados1[ob,j+1]>LSCr[j])
					{
						if(n==1)
						{
							posic[ob]_1
						}
						else
						{
							for(u in 1:m)
							{
								ranki_0
								for(ob in 1:(m*n))
								{
									if(dados1[ob,1]==u)
									{
										ranki_ranki+1
										posic[ob]_ranki
									}
								}
							} 
						}
						cat("  ",j," ",dados1[ob,1]," ",posic[ob]," ",dados1[ob,j+1])
						cat("\n")
					}
				}
			}
			cat("************************************************************************************", fill=T)
			cat("\n")
		} # fecha else
		
#################################### Fim do procedimento 5 ##################################################


################ Procedimento 6: calculo de dj's para descobrir as variaveis problematicas #################
		q2_qchisq(1-alfa,p)
		T2n_matrix(0,p,r)
		xbarras_matrix(0,p,r)
		if (problema2==1)
		{
			for(u in 1:r)	
			{
				for(j in 1:p)
				{
					# problema nesse for
					xbarras_xbarran[-j,u] # matriz (p-1)*1 
					median_media[-j] # matriz (p-1)*1 
					sn_s[-j,-j] # matriz (p-1)*(p-1)
					difn_(xbarras-median) # vetor p-1
					T2n[j,u]_n*(t(difn) %*%(solve(sn))%*% difn)
					# (1*p-1)(p-1*p-1)(p-1*1)= (1*1)
				}
			}
			#cat("T2n \n")
			#print(T2n)
			cat("\n")
			d_T2probn-T2n # matriz r*p
			qui_qchisq(1-alfa,1)
			pvalor_1-pchisq(d,1)
			pvalor_matrix(pvalor,nrow=p,byrow=F)
			cat("O numero de amostras que apresentaram problemas em T2 de Hotelling e: ",r,fill=T)
			cat("A matriz dos di's é de tal forma que as colunas mostram ")
			cat("os di's das amostras problematicas e as linhas mostram os di's das diferentes variáveis.\n")
			print(d)
			cat("Limite de corte para os di's com ",alfa," de significancia e: ",qui,fill=T)
			cat("\n A matriz com os p-valores desses di's e: \n")
			cat("\n")
			print(pvalor)
			cat("\n")
		}
		cat("****************************************************************************************", fill=T)
	}# fecha fase=2



##################################### Fim do procedimento 6 ################################################
}
 # fim do programa