Vedlegg 1 - R kode

Kapittel 2

Sette opp nødvendige pakker:

pacman::p_load(ggplot2, readxl, tidyverse, ggpubr, dplyr, hrbrthemes)

Lage fig 2.1:

prepost_eksempel_long <- as_tibble(read_excel("prepost_eksempel_long.xlsx"))

ggbarplot(prepost_eksempel_long, x = "Periode", y = "Verdi", add = c("mean"), color = "blue", fill = "lightblue")

Lage fig 2.2:

ggline(prepost_eksempel_long, x = "Periode", y = "Verdi", add = c("mean_se", "jitter"))

Lage fig 2.3:

t <- 1:24

z <- c(prepost_eksempel_long$Verdi)
  
plot(t,z, type="l", col="blue", lwd=3, xlab="Periode", ylab="Antall", xaxt="n")
axis(1, seq(0,24,2))

abline(v=12, col="red", lwd = 3)
text(15.5, 40, "Endring i prosedyre", col = "red")

Kapittel 3

Lage fig 3.1:

# Genererer tilfeldig tall for regel 1:

pacman::p_load(xlsx, ggplot2, tidyverse, ggpubr)

set.seed(91)

regel1_x <- as_tibble(rnorm(100, mean = 0, sd = 1)) %>%
  rename(x = value) %>%
  add_column(nr = 1:100) %>%
  relocate(nr)

regel1_y <- as_tibble(rnorm(100, mean = 0, sd = 1)) %>%
  rename(y = value) %>%
  add_column(nr = 1:100) %>%
  relocate(nr)

regel1 <- merge(regel1_x,regel1_y,by="nr")

regel1_plot <- ggplot(regel1, aes(x = x, y = y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 1") + xlim(-5,5) + ylim(-5,5)+ geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

# Genererer tilfeldig tall for regel 2:

set.seed(92)

regel2_x <- as_tibble(rnorm(100, mean = 0, sd = 2)) %>%
  rename(x = value) %>%
  add_column(nr = 1:100) %>%
  relocate(nr)

## Setter sd = 2 fordi variasjonen med regel 2 er dobbel av regel 1 (jfr SPC for Excel)

regel2_y <- as_tibble(rnorm(100, mean = 0, sd = 2)) %>%
  rename(y = value) %>%
  add_column(nr = 1:100) %>%
  relocate(nr)

regel2 <- merge(regel2_x,regel2_y,by="nr")

regel2_plot <- ggplot(regel2, aes(x = x, y = y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 2") + xlim(-5,5) + ylim(-5,5) + geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

# Henter tall fra Excel for regel 3:

regel3 <- read.xlsx("Funnel_blankedata_regel3.xlsx", 1)

regel3_plot <- ggplot(regel3, aes(x = Regel_3_x, y = Regel_3_y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 3") + xlim(-5,5) + ylim(-5,5)+ geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

# Henter tall fra Excel for regel 4:

regel4 <- read.xlsx("Funnel_blankedata_regel4.xlsx", 1)

regel4_plot <- ggplot(regel4, aes(x = Regel_4_x, y = Regel_4_y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 4") + xlim(-5,5) + ylim(-5,5)+ geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

# Setter plottene sammen:

ggarrange(regel1_plot, regel2_plot, regel3_plot, regel4_plot 
 + rremove("x.text"), ncol = 2, nrow = 2)

Lage fig 3.2:

regel3_2 <- read.xlsx("Funnel_blankedata_regel3_2.xlsx", 1)

regel3_2_plot <- ggplot(regel3_2, aes(x = x, y = y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 3 - annet førstetreff") + xlim(-5,5) + ylim(-5,5)+ geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

regel3_2_plot 

Lage fig 3.3:

regel4_2 <- read.xlsx("Funnel_blankedata_regel4_2.xlsx", 1)

regel4_2_plot <- ggplot(regel4_2, aes(x = x, y = y)) +
  geom_point(size = 1) + labs(x = "x", y = "y", title = "Regel 4 - annet førstetreff") + xlim(-2,10) + ylim(-10,2)+ geom_vline(xintercept = 0, col = "red") + geom_hline(yintercept = 0, col = "red")

regel4_2_plot

Lage fig 3.4:

pacman::p_load(qicharts2)

regel1x <- regel1 %>% pull(2)

regel1y <- regel1 %>% pull(3)

regel1xrun <- qic(regel1x, title = 'x-verdi ved regel 1', ylab = "verdi", xlab = "Forsøk nr.")

regel1yrun <- qic(regel1y, title = 'y-verdi ved regel 1', ylab = "verdi", xlab = "Forsøk nr.")

ggarrange(regel1xrun, regel1yrun + rremove("x.text"), ncol = 2, nrow = 1,  widths = c(1, 1))

Lage fig 3.5:

pacman::p_load(qicharts2)

regel2x <- regel2 %>% pull(2)

regel2y <- regel2 %>% pull(3)

regel2xrun <- qic(regel2x, title = 'x-verdi ved regel 2', ylab = "verdi", xlab = "Forsøk nr.")

regel2yrun <- qic(regel2y, title = 'y-verdi ved regel 2', ylab = "verdi", xlab = "Forsøk nr.")

ggarrange(regel2xrun, regel2yrun + rremove("x.text"), ncol = 2, nrow = 1,  widths = c(1, 1))

Lage fig 3.6:

pacman::p_load(qicharts2)

regel3x <- regel3 %>% pull(2)

regel3y <- regel3 %>% pull(3)

regel3xrun <- qic(regel3x, title = 'x-verdi ved regel 3', ylab = "verdi", xlab = "Forsøk nr.")

regel3yrun <- qic(regel3y, title = 'y-verdi ved regel 3', ylab = "verdi", xlab = "Forsøk nr.")

ggarrange(regel3xrun, regel3yrun + rremove("x.text"), ncol = 2, nrow = 1,  widths = c(1, 1))

Lage fig 3.7:

pacman::p_load(qicharts2)

regel4x <- regel4 %>% pull(2)

regel4y <- regel4 %>% pull(3)

regel4xrun <- qic(regel4x, title = 'x-verdi ved regel 4', ylab = "verdi", xlab = "Forsøk nr.")

regel4yrun <- qic(regel4y, title = 'y-verdi ved regel 4', ylab = "verdi", xlab = "Forsøk nr.")

ggarrange(regel4xrun, regel4yrun + rremove("x.text"), ncol = 2, nrow = 1,  widths = c(2, 2))

Lage fig 3.8:

pacman::p_load(xlsx, qicharts2, tidyverse, ggplot2, ggpubr)

toteam1 <- read.xlsx("toteam.xlsx", 1)

team1 <- qic(Team.1,
    data      = toteam1, 
    chart     = 'i',
    show.grid = TRUE,
    title     = "Team 1",
    ylab      = "Antall defekter pr uke",
    xlab      = "Uke #")

team2 <- qic(Team.2,
    data      = toteam1, 
    chart     = 'i',
    show.grid = TRUE,
    title     = "Team 2",
    ylab      = "Antall defekter pr uke",
    xlab      = "Uke #")

ggarrange(team1, team2 + rremove("x.text"), ncol = 2, nrow = 1,  widths = c(1, 1))

Lage fig 3.9:

pacman::p_load(xlsx, qicharts2, tidyverse, ggplot2, ggpubr)

toteam2 <- readxl::read_excel("toteam.xlsx") %>%
  pivot_longer(c("Team 1", "Team 2"))

qic(Nr, value,
    data = toteam2,
    facets = ~name,
    xlab = "Uke #",
    ylab = "Antall defekter pr uke",
    title = "To team - like prosesser - ulik variasjon")

Kapittel 4

Lage fig 4.1:

pacman::p_load(tidyverse)

set.seed(30)
x = rnorm(100, 179, 16)
hist(x, xlab = "Høyde", ylab = "Antall", main = "Histogram for genererte høydedata")

Lage fig 4.2:

data1 <- sample(165:175, 50, replace=TRUE)
data2 <- sample(170:180, 30, replace=TRUE)
data3 <- sample(180:185, 15, replace = TRUE)
data4 <- sample(185:190, 5, replace = TRUE)

data <- c(data1, data2, data3, data4)

hist(data, xlab = "Høyde", ylab = "Antall", main = "Histogram for genererte høydedata")

Lage fig 4.3:

pacman::p_load(ggplot2, readxl, tidyverse, ggfortify)

set.seed(100)

normalfordeling <- rnorm(100, mean = 0, sd = 1)

hist(normalfordeling, 
     main = "Genererte, normalfordelte data", 
     xlab = "x", 
     ylab = "f(x)",
     border = "black", 
     col = "gray",
     xlim = c(-4,4),
     ylim = c(0,0.5),
     las = 1, 
     probability = TRUE)

Lage fig 4.4:

set.seed(100)

normalfordeling <- rnorm(100, mean = 0, sd = 1)

hist(normalfordeling, 
     main = "Genererte, normalfordelte data", 
     xlab = "x", 
     ylab = "f(x)",
     border = "black", 
     col = "gray",
     xlim = c(-4,4),
     ylim = c(0,0.5),
     las = 1, 
     probability = TRUE)

m <- mean(normalfordeling)
std <- sqrt(var(normalfordeling))

curve(dnorm(x, mean = m, sd = std), 
      col="darkblue", lwd = 3, add = TRUE, yaxt = "n")

Lage fig 4.5:

pacman::p_load(ggplot2, readxl, tidyverse, ggfortify)

ggplot(data.frame(x = c(-4, 4)), aes(x)) +
  geom_function(fun = dnorm, colour = "darkblue", size = 1.5) +
  theme_classic() +
  scale_y_continuous(limits = c(0, 0.5), breaks = seq(0, 0.5, by = 0.1))

Lage fig 4.6:

x <- seq(-4, 4, length=200)
y <- dnorm(x)

plot(x, y, type="l", lty=1, lwd = 2, col = "red", xlab="x",
  ylab="f(x)")

x <- seq(-1,1,length=100)
y <- dnorm(x)

polygon(c(-1,x,1),c(0,y,0),col="lightblue")

abline(v=-1, col="green", lwd = 2)
text(1.3, 0.38, "1 SD", col = "black")
text(-1.35, 0.38, "-1 SD", col = "black")
abline(v=1, col="green", lwd = 2)

text(0, 0.2, "68 %", col = "black")

Lage fig 4.7:

x <- seq(-4,4,length=200)
y <- dnorm(x)

plot(x, y, type="l", lty=1, lwd = 2, col = "red", xlab="x",
  ylab="f(x)")

x <- seq(-2,2,length=200)
y <- dnorm(x)

polygon(c(-2,x,2),c(0,y,0),col="lightblue")

abline(v=-2, col="green", lwd = 2)
text(2.3, 0.38, "2 SD", col = "black")
text(-2.35, 0.38, "-2 SD", col = "black")
abline(v=2, col="green", lwd = 2)

text(0, 0.2, "95 %", col = "black")

Utregning av areal mellom to x-verdier i normalfordeling:

pnorm(2,mean=0,sd=1)-pnorm(-2,mean=0,sd=1)

Lage fig 4.8:

x <- seq(-4,4,length=200)
y <- dnorm(x)

plot(x,y,type="l",lwd=2,col="red", xlab="x",
  ylab="f(x)")

x <- seq(-3,3,length=200)
y <- dnorm(x)

polygon(c(-3,x,3),c(0,y,0),col="lightblue")

abline(v=-3, col="green", lwd = 2)
text(3.3, 0.38, "3 SD", col = "black")
text(-3.35, 0.38, "-3 SD", col = "black")
abline(v=3, col="green", lwd = 2)

text(0, 0.2, "99.7 %", col = "black")

Lage fig 4.9:

y.norm <- rnorm(n= 100000, mean = 0, sd = 1)
h <- hist(y.norm, breaks = 100, plot = F)
cuts <- cut(h$breaks, c(-Inf,-3,-2,-1,1,2,3,Inf), right = F) # right=False; sets intervals to be open on the right closed on the left
plot(h, 
     col = rep(c("white", "4","3","2","3","4", "white"))[cuts], 
     main = 'Normalfordeling', 
     xlab = '', 
     freq = F, 
     ylim = c(0,0.6))

lwd = 3
# horzintal lines
lines(x = c(2,-2), y = c(0.48,0.48), type = "l", col=3, lwd = lwd)
lines(x = c(3,-3), y = c(0.55,0.55), type = "l", col=4, lwd = lwd)
lines(x = c(1,-1), y = c(0.41,0.41), type = "l", col=2, lwd = lwd)
# vertical lines
lines(x = c(1,1), y = c(0,0.41), type = "l", col=2, lwd = lwd)
lines(x = c(-1,-1), y = c(0,0.41), type = "l", col=2, lwd = lwd)
lines(x = c(2,2), y = c(0,0.48), type = "l", col=3, lwd = lwd)
lines(x = c(-2,-2), y = c(0,0.48), type = "l", col=3, lwd = lwd)
lines(x = c(3,3), y = c(0,0.55), type = "l", col=4, lwd = lwd)
lines(x = c(-3,-3), y = c(0,0.55), type = "l", col=4, lwd = lwd)
# text
text(0, 0.44, "68%", cex = 1.5, col=2)
text(0, 0.51, "95%", cex = 1.5, col=3)
text(0, 0.58, "99.7%", cex = 1.5, col=4)

Lage fig 4.10:

success <- 0:20

plot(success,dbinom(success,size=20,prob=.5),
     type='h',
     main="Binomial distribusjon (n=20, p=0.5)",
     ylab="Sannsynlighet",
     xlab = "Suksess",
     lwd=10)

Lage fig 4.11:

success <- 0:20

plot(success,dbinom(success,size=20,prob=.2),
     type='h',
     main="Binomial distribusjon (n=20, p=0.2)",
     ylab="Sannsynlighet",
     xlab = "Suksess",
     lwd=10)

Lage fig 4.12:

set.seed(32)

terning10 <- sample(1:6, 10, replace = TRUE)

stripchart(terning10, method = "stack", offset = .5, at = 0, pch = 19,
           col = "steelblue", main = "10 terningkast", xlab = "Verdi på terning", ylab = "Antall"))

Lage fig 4.13:

set.seed(33)

terning10 <- sample(1:6, 100, replace = TRUE)

stripchart(terning10, method = "stack", offset = .5, at = 0, pch = 19,
           col = "steelblue", main = "100 terningkast", xlab = "Verdi på terning", ylab = "Antall")

Lage simulering av terningkast:

set.seed(43)
terning_runde1 <- sample(1:6, 600, replace = TRUE)
table(terning_runde1)

set.seed(44)
terning_runde2 <- sample(1:6, 600, replace = TRUE)
table(terning_runde2)

set.seed(45)
terning_runde3 <- sample(1:6, 600, replace = TRUE)
table(terning_runde3)

Lage fig 4.14:

set.seed(46)
minterning <- sample(1:6, 6000000, replace = TRUE)

options(scipen=999)
hist(minterning,
     main="Histogram for 6 000 000 terningkast",
     ylab="Antall",
     xlab = "Verdi på terning")

Lage fig 4.15:

# Grid of X-axis values
x <- 0:50

#-----------
# lambda: 5
#-----------
lambda <- 5
plot(dpois(x, lambda), type = "h", lwd = 2,
     main = "Poisson sannsynlighetsfordeling",
     ylab = "P(X = x)", xlab = "Antall hendelser")

#-----------
# lambda: 10
#-----------
lambda <- 10
lines(dpois(x, lambda), type = "h", lwd = 2, col = rgb(1,0,0, 0.7))

#-----------
# lambda: 20
#-----------

lambda <- 20
lines(dpois(x, lambda), type = "h", lwd = 2, col = rgb(0, 1, 0, 0.7))

# Legend
legend("topright", legend = c("5", "10", "20"),
       title = expression(lambda), title.adj = 0.75,
       lty = 1, col = 1:3, lwd = 2, box.lty = 0)

Lage fig 4.16:

maal <- 0:10

plot(maal, dpois(maal, lambda=2.5),
     type='h',
     main='Poissonfordeling (lambda = 2.5)',
     ylab='Sannsynlighet',
     xlab ='# Mål',
     lwd=3)

Regne ut sannsynligheter:

# Sannsynligheten for 0 mål
dpois(x = 0, lambda = 2.5)

# Sannsynligheten for 1 mål
dpois(x = 1, lambda = 2.5)
      
# Sannsynligheten for 2 mål
dpois(x = 2, lambda = 2.5)
      
# Sannsynligheten for 3 mål
dpois(x = 3, lambda = 2.5)
      
# Sannsynligheten for 4 mål
dpois(x = 4, lambda = 2.5)
      
# sannsynligheten for mellom 1 og 3 mål:
dpois(x = 1, lambda = 2.5) + dpois(x=2, lambda = 2.5) + dpois(x=3, lambda = 2.5)

Lage fig 4.17:

x_dgeom <- seq(1, 20, by = 1)

y_dgeom <- dgeom(x_dgeom, prob = 0.4) 

plot(y_dgeom,
type="l",     
main="Geometrisk fordeling for p = 0.4",
ylab="f(x)",
xlab = "x")

Lage fig 4.18:

pacman::p_load(ggpubr)

eksford <- seq(0, 20, length.out=1000)

dat1 <- data.frame(x=eksford, fx=dexp(eksford, rate=0.2)) %>%
  add_column(ID = 1:1000) %>%
  relocate(3)
             
dat2 <- data.frame(x=eksford, fx=dexp(eksford, rate=1)) %>%
  add_column(ID = 1:1000) %>%
  relocate(3)
             
dat3 <- data.frame(x=eksford, fx=dexp(eksford, rate=1.5)) %>%
  add_column(ID = 1:1000) %>%
  relocate(3)
             
dat4 <- data.frame(x=eksford, fx=dexp(eksford, rate=2)) %>%
  add_column(ID = 1:1000) %>%
  relocate(3)
             

dat1plot <- ggplot(dat1, aes(x=x, y=fx)) + geom_line() + ggtitle(expression( ~ lambda ~ " = 0.2"))
dat2plot <- ggplot(dat2, aes(x=x, y=fx)) + geom_line() + ggtitle(expression( ~ lambda ~ " = 1.0"))
dat3plot <- ggplot(dat3, aes(x=x, y=fx)) + geom_line() + ggtitle(expression( ~ lambda ~ " = 1.5"))
dat4plot <- ggplot(dat4, aes(x=x, y=fx)) + geom_line() + ggtitle(expression( ~ lambda ~ " = 2.0"))

ggarrange(dat1plot, dat2plot, dat3plot, dat4plot + rremove("x.text"), ncol = 2, nrow = 2,  widths = c(1, 1))

Kapittel 5

Lage fig 5.1:

pacman::p_load(tidyverse, ggplot2, writexl, ggpubr)

# Lage normalfordelt datasett

set.seed(89)

qqnorm <- rnorm(10000, mean=90, sd=5)

qqnorm <- as_tibble(qqnorm)

# Plotte Q-Q plott

ggqqplot(qqnorm$value) + ggtitle("Normal Q-Q plott") + labs(x = "Teoretisk forventning", y = "Data")

Lage fig 5.2:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, ggpubr)

# Lage datasett med right skew

N <- 5000
qqrightskew <- rnbinom(N, 10, .1)
 
qqrightskew <- as_tibble(qqrightskew)

# Plotte histogram og Q-Q plott

qqrighthist <- ggplot(qqrightskew, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqrightskew_plott <- ggqqplot(qqrightskew$value) + ggtitle("Normal Q-Q plott - skjevhet høyre") + labs(x = "Teoretisk forventning", y = "Data")

grid.arrange(qqrighthist, qqrightskew_plott, ncol=2)

Lage fig 5.3:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, ggpubr)

# Lage datasett med left skew

set.seed(91)

N=5000
qqleftskew <- rbeta(N,2,0.5,ncp=2)

qqleftskew <- as_tibble(qqleftskew)

# Plotte histogram og Q-Q plott"

qqlefthist <- ggplot(qqleftskew, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqlefthist <- ggplot(qqleftskew, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqleftskew_plott <- ggqqplot(qqleftskew$value) + ggtitle("Normal Q-Q plott - skjevhet venstre") + labs(x = "Teoretisk forventning", y = "Data")

grid.arrange(qqlefthist, qqleftskew_plott, ncol=2)

Lage fig 5.4:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, ggpubr)

# Lage datasett med fet hale

set.seed(14)

N=100

qqcauchy <- as_tibble(rcauchy(N, scale = 5))

# Plotte histogram og Q-Q plott

qqcauchyhist <- ggplot(qqcauchy, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqcauchyhist <- ggplot(qqcauchy, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqcauchy_plott <- ggqqplot(qqcauchy$value) + ggtitle("Normal Q-Q plott - tung hale") + labs(x = "Teoretisk forventning", y = "Data")

grid.arrange(qqcauchyhist, qqcauchy_plott, ncol=2)

Lage fig 5.5:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, ggpubr)

# Lage datasett med tynn hale

set.seed(81)

qqlt <- runif(n = 1000, min = -1, max = 1)

qqlt <- as_tibble(qqlt)

# Plotte histogram og Q-Q plott

qqlthist <- ggplot(qqlt, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqlt_plott <- ggqqplot(qqlt$value) + ggtitle("Normal Q-Q plott - lett hale") + labs(x = "Teoretisk forventning", y = "Data")

grid.arrange(qqlthist, qqlt_plott, ncol=2)

Lage fig 5.6:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, ggpubr)

# Lage bimodialt datasett

set.seed(10) 

mode1 <- rnorm(50,2,1)
mode1 <- mode1[mode1 > 0] 
mode2 <- rnorm(50,6,1)
mode2 <- mode2[mode2 > 0] 
qqbimod <- as_tibble(sort(c(mode1,mode2)))

# Plotte histogram og Q-Q plott

qqbimodhist <- ggplot(qqbimod, aes(x=value)) + geom_histogram(color="black", fill="lightblue")

qqbimod_plott <- ggplot(qqbimod, aes(sample = value)) +
  stat_qq() +
  stat_qq_line() +
  ggtitle(" Normal Q-Q plott - bimodial") + labs(x = "Teoretisk forventning", y = "Data")

grid.arrange(qqbimodhist, qqbimod_plott, ncol=2)

Kjøre Anderson-Darling test for normalitet:

pacman::p_load(nortest, readxl, tidyverse)

addata <- as_tibble(read_excel("Anderson-Darling_raw.xlsx"))

ad.test(addata$Values)

Lage fig 5.7:

pacman::p_load(gridExtra, writexl, ggplot2, tidyverse, readxl)

addata2 <- as_tibble(read_excel("Anderson-Darling_raw.xlsx"))

ggplot(addata2, aes(sample = Values)) +
  stat_qq() +
  stat_qq_line() +
  ggtitle(" Normal Q-Q plott - A-D data") + labs(x = "Teoretisk forventning", y = "Data")

Kjøre Anderson-Darling test for normalitet på normalfordelte data:

pacman::p_load(nortest, readxl, tidyverse)

addata3 <- as_tibble(read_excel("QQ_norm.xlsx"))

ad.test(addata3$value)

Kjøre statistiske tester:

pacman::p_load(nortest, readxl, tidyverse, tseries)

options(scipen=999)

addata5 <- as_tibble(read_excel("Anderson-Darling_raw.xlsx"))

ks.test(addata5, "pnorm")

shapiro.test(addata5$Values)

cvm.test(addata$Values)

Jarque-Bera test:

addata6 <- as_tibble(read_excel("Anderson-Darling_raw.xlsx"))

pacman::p_load(tseries, normtest)

jarque.bera.test(addata6$Values)

ajb.norm.test(addata6$Values, nrepl=2000)

Lage fig 5.8:

pacman::p_load(ggpubr, tidyverse)

addata7 <- as_tibble(read_excel("Anderson-Darling_raw.xlsx"))

ggqqplot(addata7$Values)

Kapittel 6

Lage fig 6.1:

pacman::p_load(qicharts2, ggplot2)

set.seed(21)

eksempelrun <- rnorm(24, 16)

serie1 <- qic(eksempelrun, title = "Eksempel på seriediagram - normalfordelte genererte tall", ylab = "Verdi", xlab = "Hendelse",   method = "anhoej")

serie1 + geom_point(size = 2)

Lage fig 6.2:

pacman::p_load(qicharts2, ggplot2)

set.seed(43)

eksempelrun[13:24] <- rpois(12, 24)  

serie2 <- qic(eksempelrun, title = "Eksempel på seriediagram - modifiserte genererte tall", ylab = "Verdi", xlab = "Hendelse", method = "anhoej")

serie2 + geom_point(size = 2)

Lage fig 6.3:

pacman::p_load(qicharts2, ggplot2)

serie3 <- qic(eksempelrun, title = "Eksempel på seriediagram - endring i prosess", ylab = "Verdi", xlab = "Hendelse", method = "anhoej", part = 14)

serie3 + geom_point(size = 2)

Kapittel 7

Lage fig 7.1:

pacman::p_load(qicharts2, tidyverse)

set.seed(81)

kontr_eks <- (rnorm(24))

kontr1 <- qic(kontr_eks, chart = 'i', title = "Eksempel på kontrolldiagram", subtitle = "Tilfeldig genererte tall", ylab = "Verdi", xlab = "Hendelse") 

kontr1 + geom_point(size = 2)

Lage fig 7.2:

pacman::p_load(qcc, tidyverse)

set.seed(81)

kontr_eks <- (rnorm(24))

kontr_eks <- as_tibble(kontr_eks)

kontr2 <- qcc(kontr_eks, type = "xbar.one", nsigmas = 3)

plot(kontr2, title = "Eksempel på kontrolldiagram - Tilfeldig genererte tall", ylab = "Verdi", xlab = "Hendelse")

Lage fig 7.3:

pacman::p_load(qcc, tidyverse)

set.seed(64)

mode1x <- rnorm(15,2,1)
mode1x <- mode1x[mode1x > 0] 
mode2x <- rnorm(15,3,2)
mode2x <- mode2x[mode2x > 0] 
kontr_eks2 <- (sort(c(mode1x,mode2x)))

kontr_eks2 <- as_tibble(kontr_eks2)

qq <- qcc(kontr_eks2, type = "xbar.one", nsigmas = 3)

plot(qq, title = "Eksempel på kontrolldiagram - Tilfeldig genererte tall", ylab = "Verdi", xlab = "Hendelse")

Lage fig 7.4:

pacman::p_load(qcc, tidyverse, readxl)

pdiagr <- as_tibble(read_excel("P_chart.xlsx"))

p_chart <- with(pdiagr, qcc(pdiagr$Keisersnitt, pdiagr$Antall_fødsler, type = "p"))

plot(p_chart, title = "p-diagram: Andel keisersnitt", xlab = "Måned", ylab = "Andel")

Lage fig 7.5:

pacman::p_load(qcc, tidyverse, readxl)

Lpdiagr <- as_tibble(read_excel("Laneyp.xlsx"))

Lp_chart <- with(Lpdiagr, qcc(Lpdiagr$Pr_telefon, Lpdiagr$Medlemmer, type = "p"))

plot(Lp_chart,  title = "p-diagram: Andel kommunisert pr telefon", xlab = "Måned", ylab = "Andel")

Lage fig 7.6:

pacman::p_load(qcc, tidyverse, readxl)

npdiagr <- as_tibble(read_excel("np_diagram.xlsx"))

np_chart <- with(npdiagr, qcc(npdiagr$Feil, npdiagr$n, type = "np"))

plot(np_chart, title = "np-diagram: Andel feil", xlab = "Uke", ylab = "Andel")

Lage fig 7.7:

pacman::p_load(qcc, tidyverse, readxl)

udiagr <- as_tibble(read_excel("u_diagram.xlsx"))

u_chart <- with(udiagr, qcc(udiagr$Pasientfall, udiagr$Pasientdager, type = "u"))

plot(u_chart,  title = "u-diagram: Andel feil", xlab = "Uke", ylab = "Pasientfall")

Lage fig 7.8:

pacman::p_load(qcc, tidyverse, readxl)

cdiagr <- as_tibble(read_excel("cdiagram.xlsx"))

c_chart <- with(cdiagr, qcc(cdiagr$Feilmedisinering, type = "c"))

plot(c_chart, title = "c-diagram", xlab = "Måned", ylab = "Feilmedisinering")

Lage fig 7.9:

pacman::p_load(qcc, tidyverse, readxl)

imrdiagr <- as_tibble(read_excel("imr_diagram.xlsx"))

imr_chart <- with(imrdiagr, qcc(imrdiagr$Fødselsvekt, type = "xbar.one"))

plot(imr_chart, title = "IMR-diagram", xlab = "Baby nr", ylab = "Fødselsvekt")

Lage fig 7.10:

pacman::p_load(qicharts2, tidyverse)

set.seed(43)

datoer  <- seq(as.Date('2020-1-1'), as.Date('2020-12-31'), by = 'day')
hendelser <- sort(sample(datoer, 24))

t_diagram_data <- c(NA, diff(hendelser))

t_diagram <- qic(t_diagram_data,
    chart = "t",
    title = "T-diagram",
    xlab = "Hendelse nr.",
    ylab = "Dager")

t_diagram + geom_point()

Vedlegg 2

pacman::p_load(kableExtra)

n <- 10:100
hendelser <- data.frame(
  "Antall observasjoner"                = n,
  "Øvre grense for serie"               = round(log2(n) + 3),
  "Nedre grense for antall krysninger"  = qbinom(0.05, n - 1, 0.5),
  check.names = FALSE)

kbl(hendelser, booktabs = T,  longtable = T, caption = "Kritiske verdier", align = 'c') %>%
  kable_styling(latex_options = "striped") %>%
  kable_styling(latex_options = "repeat_header") %>%
  row_spec(0, bold = T)

Vedlegg 3

Lage generisk grafisk framstilling av Chebyshevs teorem

pacman::p_load(tidyverse)

# For å lage eksempelet lager vi to normalfordelte datasett med ulik gjennomsnitt og standardavvik som vi slår sammen

set.seed(10) 

mode1 <- rnorm(1000,2,1)
mode1 <- mode1[mode1 > 0] 
mode2 <- rnorm(1000,6,2)
mode2 <- mode2[mode2 > 0] 
modex2 <- as_tibble(sort(c(mode1,mode2)))

# Deretter setter vi dataene inn i et diagram

ggplot(modex2, aes(x=value)) + 
  geom_density() +
  ggtitle("Eksempel: Bimodal distribusjon", subtitle = "Eksemplet er kun illustrativt, ikke nøyaktig eller basert på reelle data") +
  labs(x = "", y = "") +
  theme_classic() +
  theme(legend.position = "none") + 
  
  scale_x_discrete(labels = NULL, breaks = NULL) +
  labs(x = "") +
  xlim(-0.02, 8) +
  ylim(-0.02,0.3) +
 
  annotate("segment", x = 4, y = 0, xend = 4, yend = 0.15, linetype = "dashed", color = "red") + 
  annotate('text', x = 4, y = -0.015, label = "bar(x)", parse = TRUE, size = 5) +
  
  annotate("segment", x = 7.8, y = 0, xend = 7.8, yend = 0.275, color = "darkgreen") + 
  annotate('text', x = 7.8, y = -0.015, label = "bar(x) ~ + ~ 3 ~ sd", parse = TRUE, size = 5) +
  
  annotate("segment", x = 0.2, y = 0, xend = 0.2, yend = 0.275, color = "darkgreen") + 
  annotate('text', x = 0.3, y = -0.015, label = "bar(x)~-~3~sd", parse = TRUE, size = 5) +
  
  annotate("segment", x = 2.1, y = 0, xend = 2.1, yend = 0.235, color = "blue") + 
  annotate('text', x = 2.1, y = -0.015, label = "bar(x)~-~2~sd", parse = TRUE, size = 5) +
  
  annotate("segment", x = 5.9, y = 0, xend = 5.9, yend = 0.235, color = "blue") + 
  annotate('text', x = 5.9, y = -0.015, label = "bar(x)~+~2~sd", parse = TRUE, size = 5) +
  
  annotate("text", x=4, y=0.23, label="Minst 75 %", color = "blue") +
  annotate("segment", x = 3.35, y = 0.23, xend = 2.25, yend = 0.23, color = "blue") + 
  annotate("segment", x = 4.6, y = 0.23, xend = 5.8, yend = 0.23, color = "blue") +
  
  annotate("text", x=4, y=0.27, label="Minst 88.89 %", color = "darkgreen") +
  annotate("segment",x = 0.3, y = 0.27, xend = 3.2, yend = 0.27, color = "darkgreen") +
  annotate("segment", x = 4.8, y = 0.27, xend = 7.7, yend = 0.27, color = "darkgreen")

Lage tabell over k og prosent

pacman::p_load(tidyverse, kableExtra, knitr)
k <- seq(1,4,by = 0.1)
auc <- 1-(1/k^2)
auc.percent <- round(auc*100)
cheb_table <- as_tibble(cbind(k,auc.percent))

kbl(cheb_table) %>%
  kable_styling(bootstrap_options = "striped", full_width = F) %>%
  kable_paper() %>%
  scroll_box(height = "200px")

Lage graf over k og prosent

plot(k, 
     auc.percent, 
     col = 'blue', 
     pch = 10, 
     xlab = 'k', 
     ylab = 'Prosent', 
     main = 'Chebyshevs teorem' )

abline(v=2, col="red", lwd = 1)
abline(h=75, col="red", lwd = 1)

Vedlegg 4

Lage ikke-normalfordelte data og plotte Q-Q plott

pacman::p_load(ggplot2, tidyverse, magrittr, dplyr, truncnorm)

nn <- 1e4
set.seed(1)
sims <- as_tibble(c(rtruncnorm(nn/2, a=1, b=5, mean=2, sd=.5),
                    rtruncnorm(nn/2, a=1, b=5, mean=4, sd=.5)))

ggplot(sims, aes(sample=value)) + stat_qq() + stat_qq_line(col = "red")

Lage første histogram (populasjonen)

pacman::p_load(ggplot2, tidyverse, magrittr, dplyr, truncnorm)

nn <- 1e4
set.seed(1)
sims <- c(rtruncnorm(nn/2, a=1, b=5, mean=2, sd=.5),
                    rtruncnorm(nn/2, a=1, b=5, mean=4, sd=.5))

mySD <- as.character(abs(as.integer((sims - mean(sims)) / sd(sims))))
myDF <- data.frame(sims, mySD)
xAxis <- as.integer(max(abs(sims)))
mu <- round(mean(sims),2)
sd <- round(sd(sims),2)
myBin <- sd/10
ggplot(myDF, aes(sims)) +
  geom_histogram(aes(fill = mySD), binwidth = myBin, col="black", size=.1) +  # change binwidth
  labs(x="x", y="Frequency") + 
  labs(title="Histogram av generert bimodal distribusjon",subtitle=paste0(  "Gj.snitt = ", mu, ", sd = ", sd)) +
  scale_x_continuous(breaks = seq(mu-sd*5, mu+sd*5, sd)) +
  theme_bw() +
  theme(plot.title = element_text(hjust = 0.5)) +
  guides(fill=guide_legend(expression(sigma))) 

Lage histogram over 100 utvalg av 30

n <- 100
sampSize <- 30
xbar <- rep(NA, n)
for (i in 1:n) {
  mysamp <- sample(sims, size = sampSize)
  xbar[i] <- mean(mysamp)
}

mySD<-as.character( abs(as.integer((xbar - mean(xbar)) / sd(xbar) )))
myDF<-data.frame(xbar,mySD)
xAxis<-as.integer(max(abs(xbar)))
mu<-round(mean(xbar),2)
sd<-round(sd(xbar),2)
myBin<-sd/10
ggplot(myDF, aes(xbar)) +
  geom_histogram(aes(fill = mySD), binwidth = myBin, col="black", size=.1) +  # change binwidth
  labs(x="x", y="Frequency") + 
  labs(title="Histogram for genererte utvalg",
  subtitle=paste0(  "mean = ", mu, ", sd = ", sd, ", Utvalgsstørrelse = ",sampSize,", Antall utvalg = ",n))+
  scale_x_continuous(breaks = seq(mu-sd*5, mu+sd*5, sd))+
  theme_bw()+
  theme(plot.title = element_text(hjust = 0.5))+
  guides(fill=guide_legend(expression(sigma)))+
  geom_density(aes(y=..count../90))

Lage histogram over 1000 utvalg av 30

n<-1000
sampSize<-30
xbar <- rep(NA, n)
for (i in 1:n) {
  mysamp <- sample(sims, size = sampSize)
  xbar[i] <- mean(mysamp)
}

mySD<-as.character( abs(as.integer((xbar - mean(xbar)) / sd(xbar) )))
myDF<-data.frame(xbar,mySD)
xAxis<-as.integer(max(abs(xbar)))
mu<-round(mean(xbar),2)
sd<-round(sd(xbar),2)
myBin<-sd/10
ggplot(myDF, aes(xbar)) +
  geom_histogram(aes(fill = mySD), binwidth = myBin, col="black", size=.1) +  # change binwidth
  labs(x="x", y="Frequency") + 
  labs(title="Histogram for genererte utvalg",
  subtitle=paste0(  "Gj.snitt = ", mu, ", sd = ", sd, ", Utvalgsstørrelse = ",sampSize,", Antall utvalg = ",n))+
  scale_x_continuous(breaks = seq(mu-sd*5, mu+sd*5, sd))+
  theme_bw()+
  theme(plot.title = element_text(hjust = 0.5))+
  guides(fill=guide_legend(expression(sigma)))+
  geom_density(aes(y=..count../90))

Lage histogram over 1000 utvalg av 30

n<-10000
sampSize<-30
xbar <- rep(NA, n)
for (i in 1:n) {
  mysamp <- sample(sims, size = sampSize)
  xbar[i] <- mean(mysamp)
}

mySD<-as.character( abs(as.integer((xbar - mean(xbar)) / sd(xbar) )))
myDF<-data.frame(xbar,mySD)
xAxis<-as.integer(max(abs(xbar)))
mu<-round(mean(xbar),2)
sd<-round(sd(xbar),2)
myBin<-sd/10
ggplot(myDF, aes(xbar)) +
  geom_histogram(aes(fill = mySD), binwidth = myBin, col="black", size=.1) +  # change binwidth
  labs(x="x", y="Frequency") + 
  labs(title="Histogram for genererte utvalg",
  subtitle=paste0(  "Gj.snitt = ", mu, ", sd = ", sd, ", Utvalgsstørrelsee = ",sampSize,", Antall utvalg = ",n))+
  scale_x_continuous(breaks = seq(mu-sd*5, mu+sd*5, sd))+
  theme_bw()+
  theme(plot.title = element_text(hjust = 0.5))+
  guides(fill=guide_legend(expression(sigma)))+
  geom_density(aes(y=..count../90))

Vedlegg 5

pacman::p_load(flextable, SixSigma, officer, magrittr)

nmax <- 25
n <- 2:nmax

d2 <- sapply(2:nmax, ss.cc.getd2)
d3 <- sapply(2:nmax, ss.cc.getd3)
c4 <- sapply(2:nmax, ss.cc.getc4)
A2 <- 3/(d2*sqrt(n))
D3 <- sapply(1:(nmax-1), function(x){
max(c(0, 1 - 3*(d3[x]/d2[x])))})
D4 <- (1 + 3*(d3/d2))
B3 <- sapply(1:(nmax-1), function(x){
max(0, 1 - 3*(sqrt(1-c4[x]^2)/c4[x]))})
B4 <- 1 + 3*(sqrt(1-c4^2)/c4)

constdf <- data.frame(n, A2, d2, d3, c4, 
D3, D4, B3, B4)

set_flextable_defaults(big.mark = " ", 
  font.size = 10, 
  theme_fun = theme_vanilla,
  padding.bottom = 6, 
  padding.top = 6,
  padding.left = 6,
  padding.right = 6,
  background.color = "#EFEFEF")

constdf1 <- flextable(constdf)

constdf1 <- align(constdf1, align = "center", part = "all")

constdf1 <- add_header_lines(constdf1, values = "Tabell med konstanter for kontrolldiagram")

constdf1