Nothing
# This file exists to generation the simulation tabulated results for detection
# times and the histograms for false alarm rates
###### Multi-State Simulation ###############################################
# Read in the latest simulation data
library(xts)
library(readr)
library(reshape2)
library(tidyverse)
detection_times2 <- read_csv("~/Users/gabriel_odom/Box Sync/Consulting/Dr. Hering/MV_Process_Control/MVSPC/Simulation_Data/detection_times_ms_1000_20170323.csv")
detection_times_ms <- read_csv("C:/Users/gabriel_odom/Box Sync/Consulting/Dr. Hering/MV_Process_Control/MVSPC/Simulation_Data/detection_times_ms_1000_20170323.csv")
# Table out summary of detection times
mspSummary <- function(vec,
quants = c(0, 0.05, 0.25, 0.5, 0.75, 0.95, 1),
oobCode = Inf){
# Take in a vector, find the proportion out of bounds (NA, NaN, or Inf), and
# calculate quantile statistics (min, 0.05, 0.25, 0.5, mean, 0.75, 0.95, max)
# and spread statistics (IQR, standard deviation). Return this info as a named
# vector
vec <- vec %>% as.matrix() %>% as.vector()
### Out of Bounds Observations ###
# subset out the out of bounds observations
vecClean <- vec[!(vec %in% oobCode)]
oobProp <- 1 - length(vecClean) / length(vec)
### Summary Statistics ###
# Make sure the quantiles necessary for the IQR are included
quantsFull <- union(quants, c(0.25, 0.75))
# Calculate summaries
vecQuants <- quantile(vecClean, probs = quantsFull)
vecIQR <- vecQuants["75%"] - vecQuants["25%"]
names(vecIQR) <- NULL
vecMean <- mean(vecClean)
vecSD <- sd(vecClean)
### Return ###
returnVec <- c(vecQuants,
"IQR" = vecIQR,
"Mean" = vecMean,
"StdDev" = vecSD,
"Prop_OOB" = oobProp)
}
# Fault 1A
ms_detectionTimes5 %>%
filter(Faults == "A1") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 1B
ms_detectionTimes5 %>%
filter(Faults == "B1") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 1C (State 3)
# 8500 / 60 < 142, so the next state start at 142 * 60 = 8520. The fault state
# is one after that, so another 60 minutes. This fault starts at 80, so the
# earliest we should detect a fault is at 80. We need to remove all "alarms"
# that happen before 80, because these are false.
ms_detectionTimes4 %>%
filter(Faults == "C1") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 2A
ms_detectionTimes5 %>%
filter(Faults == "A2") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 2B
ms_detectionTimes5 %>%
filter(Faults == "B2") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 2C (State 2)
# 8500 / 60 < 142, so the next state start at 142 * 60 = 8520. This fault starts
# at 20, so the earliest we should detect a fault is at 20.
ms_detectionTimes4 %>%
filter(Faults == "C2") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 3A
ms_detectionTimes5 %>%
filter(Faults == "A3") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 3B
ms_detectionTimes5 %>%
filter(Faults == "B3") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 3C (State 2)
# 8500 / 60 < 142, so the next state start at 142 * 60 = 8520. This fault starts
# at 20, so the earliest we should detect a fault is at 20.
ms_detectionTimes4 %>%
filter(Faults == "C3") %>%
select(-Faults) %>%
sapply(mspSummary)
### False Alarm Rates ###
FA_rates2 <- read_csv("~/Box Sync/Consulting/Dr. Hering/MV_Process_Control/MVSPC/Simulation_Data/false_alarm_rates_1000_20170314.csv")
ms_FA_rates5$MSAD_SPE %>%
hist(xlim = c(0, 0.06), main = "MSAD SPE False Alarm Rates")
ms_FA_rates5$AD_T2 %>%
hist(xlim = c(0, 0.06), main = "AD T2 False Alarm Rates")
###### Single-State Simulation ##############################################
### Detection Times ###
# Read Detection Times Data
# Fault 1A
ss_detectionTimes %>%
filter(Faults == "A1") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 1B
ss_detectionTimes %>%
filter(Faults == "B1") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 2A
ss_detectionTimes %>%
filter(Faults == "A2") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 2B
ss_detectionTimes %>%
filter(Faults == "B2") %>%
select(-Faults) %>%
sapply(mspSummary)
# Fault 3A
ss_detectionTimes %>%
filter(Faults == "A3") %>%
select(-Faults) %>%
sapply(mspSummary) %>% round(3)
# Fault 3B
ss_detectionTimes %>%
filter(Faults == "B3") %>%
select(-Faults) %>%
sapply(mspSummary)
### False Alarm Rates ###
# Read in False Alarm Rates Data
ss_FA_rates[,1] %>% hist(xlim = c(0, 0.03), main = "MSAD SPE False Alarm Rates")
ss_FA_rates[,3] %>% hist(xlim = c(0, 0.03), main = "AD SPE False Alarm Rates")
ss_FA_rates[,2] %>% hist(xlim = c(0.05, 0.2), main = "MSAD T2 False Alarm Rates")
ss_FA_rates[,4] %>% hist(xlim = c(0.05, 0.2), main = "AD T2 False Alarm Rates")
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.