# Sete Ferramentas Básicas
#
# 1 -  Diagrama de Ishikawa
# 2 -  Lista de verificação
# 3 -  Carta de controle
# 4 -  Histograma
# 5 -  Diagrama de Pareto
# 6 -  Diagrama de dispersão
# 7 -  Fluxograma
# 7'-  Estratificação
#
# Pacotes: qcc, SixSigma, DiagrammeR
#

# Walter Ponge-Ferreira
# Cotia, 1º de Abril de 2018

library(qcc)
library(SixSigma)
library(moments)

setwd("~/Documents/Data/Pasta Técnica/R/PME3463 Aula 03")
rm(list=ls())

# Diagrama de Ishikawa

cManpower <- c("Instalação","Manutenção","Usuário","Limpeza")

cMaterials <- c("Energia elétrica","Lâmpadas","Reatores","Interruptores")

cMachines <- c("Relede Elétrica","Leiaute da Iluninação", "Janelas e iluminação natural")

cMethods <- c("Aulas com lousa","Projeções","Aulas de laboratório","Manuseio dos equipamento")

cMeasurements <- c("Nível de iluminamento", "Consumo de Energia Elétrica", "Relógio", "Calendário", "Bússula", "GPS")

cEnvironment <- c("Localização", "Clima", "Diurno/Noturno", "Estação do Ano")

cGroups <- c("Mão de Obra", "Materiais", "Máquina", "Método", "Metrologia", "Meio Ambiente")

cEffect <- "Iluminação em Sala de Aula"

cause.and.effect(
  cause = list(Mão_de_Obra = cManpower,
               Materiais = cMaterials,
               Máquinas = cMachines,
               Método = cMethods,
               Metrologia = cMeasurements,
               Meio_Ambinete = cEnvironment) ,
  effect = cEffect)


ss.ceDiag(
  effect = cEffect,
  causes.gr = cGroups,
  causes = list(cManpower, cMaterials, cMachines, cMethods, cMeasurements, cEnvironment),
  main = "Diagrama de Causa-e-Efeito",
  sub = "Iluminação em Sala de Aula")

# Lista de Verificação

data_checkSheet <- rbind(
  data.frame(Group = "Mão de Obra",
             Cause = cManpower),
  data.frame(Group = "Máquinas",
             Cause = cMachines),
  data.frame(Group = "Materiais",
             Cause = cMaterials),
  data.frame(Group = "Métodos",
             Cause = cMethods),
  data.frame(Group = "Metrologia",
             Cause = cMeasurements),
  data.frame(Group = "Meio Ambiente",
             Cause = cEnvironment)
)

data_checkSheet$A_supplier <- NA
data_checkSheet$B_supplier <- NA
data_checkSheet$C_supplier <- NA

str(data_checkSheet)
head(data_checkSheet)

write.csv2(data_checkSheet, file = "CheckSheet.csv")

# Carta de Controle

pdensity <- c(10.6817, 10.6040, 10.5709, 10.7858,
              10.7668, 10.8101, 10.6905, 10.6079,
              10.5724, 10.7736, 11.0921, 11.1023,
              11.0934, 10.8530, 10.6774, 10.6712,
              10.6935, 10.5669, 10.8002, 10.7607,
              10.5470, 10.5555, 10.5705, 10.7723)

myControlChart <- qcc(data = pdensity,
                      type = "xbar.one")

summary(myControlChart)

myControlChart$violations

# Histograma

hist(pdensity)

opar <- par(no.readonly = TRUE)
par(bg = "gray95")

hist(pdensity,
     main = "Histogram of pellets density - Samplo #25",
     sub = "Data from ceramic process",
     xlab = expression("Density (g"/"cm"^3*")"),
     col = "steelblue",
     border = "white",
     lwd = 2,
     las = 1,
     bg = "gray")
par(opar)

library(lattice)
histogram(pdensity,
          xlab = expression("Pellets density(g"/"cm"^3*")"),
          ylab = "Probability density",
          type = "density",
          panel = function(x, ...){
            panel.histogram(x, ...)
            panel.mathdensity(dmath = dnorm,
                              col = "black",
                              lwd = 3,
                              args = list(mean = mean(x),
                                          sd = sd(x)))
          })

library(ggplot2)
ggplot(data = data.frame(pdensity),
       aes(x = pdensity)) +
  geom_histogram(fill = "seagreen",
                 colour = "lightgoldenrodyellow",
                 binwidth = 0.2) +
  labs(title = "Histogram",
       x = expression("Density ("*g/cm^3*")"),
       y = "Frequency")

# Diagrama de Pareto

data_checkSheet$A_supplier <- c(2, 0, 0, 2, 1,
                                7, 1, 3, 6, 0,
                                1, 2, 0, 2, 1,
                                0, 0, 3, 1, 2,
                                1, 1, 0, 2, 4)
data_checkSheet$B_supplier <- c(0, 0, 1, 1, 2,
                                1, 12, 1, 1, 1,
                                0, 0, 3, 0, 1,
                                1, 2, 1, 0, 0,
                                1, 0, 4, 0, 1)
data_checkSheet$C_supplier <- c(0, 1, 0, 6, 0,
                                2, 2, 4, 3, 0,
                                1, 0, 2, 0, 0,
                                1, 1, 1, 1, 1,
                                0, 3, 4, 0, 0)
data_checkSheet$Total <- data_checkSheet$A_supplier +
  data_checkSheet$B_supplier +
  data_checkSheet$C_supplier

data_checkSheet

barplot(heigh = data_checkSheet$Total,
        names.arg = data_checkSheet$Cause)

data_pareto <- data_checkSheet[order(
  data_checkSheet$Total,
  decreasing = TRUE), ]

opar <- par(no.readonly = TRUE)
par(mar = c(8, 4, 4, 2) + 0.1)
barplot(height = data_pareto$Total,
        names.arg = data_pareto$Cause,
        las = 2,
        main = "Parreto Chart for total cases")
par(opar)

library(qcc)
data_pareto2 <- data_pareto$Total
names(data_pareto2) <- data_pareto$Cause
pareto.chart(data_pareto2,
             main = "Out-of-control causes")

library(qualityTools)
paretoChart(x = data_pareto2,
            main = "Out-of-control causes")

library(qicharts)
spreadvector <- rep(names(data_pareto2),
                    times = data_pareto2)
paretochart(spreadvector)


# Exemplo com dados dos aneis metálicos

rm(list=ls())

dados <- read.csv2("Amostra_02_Formulario_02.csv")

nomes_s <- c("oval","achatado","bicudo","disforme","assimetrica","curva","chamfrada","longa.curta")
s <- dados[,nomes_s]

ss <- apply(s,2,sum)

barplot(ss, main = "Lote Piloto", ylab = "Nº de defeitos")

ss <- sort(ss, decreasing = TRUE)
sc <- cumsum(ss)
barplot(ss, main = "Lote Piloto", ylab = "Nº de defeitos")
barplot(ss/sum(ss)*100, main = "Lote Piloto", ylab = "proporção (%)")

pareto.chart(ss, main = "Lote Piloto", ylab = "frequência", ylab2 = "proporção acumulada")

# Exemplo Wikipedia

rm(list=ls())
s <- c(82.70, 11.75, 2.30, 1.85, 1.40)
names(s) <- c("Richest 20%","Second 20%", "Third 20%", "Fourth 20%", "Poorest 20%")
s

pareto.chart(s, main = "Distribuição do PIB Mundia (1989)", ylab = "frequência", ylab2 = "proporção acumulada")

# Histograma

rm(list=ls())

dados <- read.csv2("Amostra_02_Formulario_01.csv")
str(dados)

nomes_d <- c("d01","d02","d03","d04","d05","d06","d07","d08","d09","d10")
d <- dados[nomes_d]
pc <- as.factor(seq(1,length(d[,1]),1))
d <- cbind(pc, d)
str(d)
dm <- apply(d[,nomes_d],1,mean)

library(reshape2)
dl <- melt(d, id = "pc")
hist(dl[dl$variable=="d01","value"])

attach(dl)
n <- length(value)
k <- ceiling(log2(n))
R <- diff(range(value))
D <- signif(R/k,1)
bins <- seq(trunc(min(value)),ceiling(max(value)),by = D)
detach(dl)
hist(dl$value, breaks = bins, main = "Histograma",
     sub = "Diâmetro interno dos aneis",
     xlab = "diâmetro - d / mm",
     ylab = "frequência - f")
text(min(bins)+1/8*diff(range(bins)),20,paste("n =",length(dl$value), sep = " "))



nomes_m <- c("m01","m02","m03","m04","m05","m06","m07","m08","m09","m10")
m <- dados[nomes_m]
pc <- as.factor(seq(1,length(d[,1]),1))
m <- cbind(pc, m)
str(d)

ml <- melt(m, id = "pc")
hist(ml[ml$variable=="m01","value"])
str(ml)
mm <- apply(m[,nomes_m],1,mean)


# Diagrama de Dispersão

rm(list=ls())

dados <- read.csv2("Amostra_02_Formulario_01.csv")
str(dados)

nomes_d <- c("d01","d02","d03","d04","d05","d06","d07","d08","d09","d10")
d <- dados[nomes_d]
pc <- as.factor(seq(1,length(d[,1]),1))
d <- cbind(pc, d)
str(d)
dm <- apply(d[,nomes_d],1,mean)

library(reshape2)
dl <- melt(d, id = "pc")

nomes_m <- c("m01","m02","m03","m04","m05","m06","m07","m08","m09","m10")
m <- dados[nomes_m]
pc <- as.factor(seq(1,length(d[,1]),1))
m <- cbind(pc, m)
str(d)

ml <- melt(m, id = "pc")
mm <- apply(m[,nomes_m],1,mean)

library(gplots)
plotmeans(as.matrix(d[,nomes_d]) ~ d[,1], pch = 16, main = "Lote Piloto", xlab = "Peça", ylab = "diâmetro - d / mm")

boxplot(t(d[,nomes_d]), main = "Lote Piloto", xlab = "Peça", ylab = "diâmetro - d / mm")

plotmeans(as.matrix(m[,nomes_m]) ~ m[,1], pch = 16, main = "Lote Piloto", xlab = "Peça", ylab = "massa - m / g")

boxplot(t(m[,nomes_m]), main = "Lote Piloto", xlab = "Peça", ylab = "massa - m / g")

plot(mm ~ dm, pch = 16, main = "Lote Piloto", xlab = "diâmetro - d /mm", ylab = "massa - m / g")
abline(h=mean(mm), lty=3, col="red")
abline(v=mean(dm), lty=3, col="red")
text(19.9,0.718,"2 pontos")
text(20.7,0.718,"4 pontos")
text(19.9,0.698,"3 pontos")
text(20.7,0.698,"1 ponto")
mod <- lm(mm ~ dm)
abline(mod, lty=5, col="blue")
cor(mm,dm)


# Estratificação

rm(list=ls())

load("ASQ Stratification")

str(dados)

d1 <- dados[,1:2]
d2 <- dados[,3:4]
d3 <- dados[,5:6]

names(d1) <- c("X","Y")
names(d2) <- c("X","Y")
names(d3) <- c("X","Y")

d1$id <- "Reator 1"
d2$id <- "Reator 2"
d3$id <- "Reator 3"

d <- rbind(d1,d2,d3)
d$id <- as.factor(d$id)
d <- na.omit(d)

plot(d$X, d$Y, pch = 16, cex = 1.5, col = "blue",
     main = "Contaminação do Aço - Z44",
     sub = "(sem estraficação)",
     ylab = "pureza (%)",
     xlab = "contaminação (ppm)")

boxplot(d$X, main = "Produção de Aço - Z44", sub = "(sem estratificação)", ylab = "contaminação (ppm)")
boxplot(d$Y, main = "Produção de Aço - Z44", sub = "(sem estratificação)", ylab = "pureza (%)")

boxplot(d$X ~ d$id, main = "Produção de Aço - Z44", sub = "(com estratificação)", ylab = "contaminação (ppm)")
boxplot(d$Y ~ d$id, main = "Produção de Aço - Z44", sub = "(com estratificação)", ylab = "pureza (%)")

plot(d$X, d$Y, col = d$id, pch = 16, cex = 1.5,
     main = "Contaminação do Aço - Z44",
     sub = "(com estratificação)",
     ylab = "pureza (%)",
     xlab = "contaminação (ppm)")
legend("topright",as.character(sort(unique(d$id))), fill = sort(unique(d$id)))