Kadrovska služba, ki skrbi za izobraževanje zaposlenih, je poskušala določiti optimalno število izobraževanj za skupino zaposlenih, ki delajo na proizvodni liniji. Naključno so izbrali 38 novih zaposlenih in opazovali število izdelkov, ki so jih proizvedli na mesec. Nato so jih poslali na prvo usposabljanje in ponovili opazovanje števila izdelanih izdelkov. Poslali so jih tudi na drugo in tretje usposabljanje. Koliko izobraževanj naj podjetje ponudi zaposlenim na podlagi teh podatkov?

ENOTA: zaposlen. n= 38, spremenljivka = število izdelkov po usposabljanju

library(readxl)
podatki <- read_xlsx("./Produktivnost.xlsx")
podatki <- as.data.frame(podatki)
head(podatki)
##   ID Zacetek PU1 PU2 PU3
## 1  1      85 133 136 134
## 2  2      90 123 138 135
## 3  3      96 120 135 141
## 4  4      93 129 132 134
## 5  5      84 117 138 131
## 6  6      91 128 132 141

Imamo širok format podatkov.

Opis:

Analiza variance za odvisne vzorce: rANOVA. PREDPOSTAVKE: 1. Število izdelkov mora biti porazdeljeno normalno na začetku, po 1. 2. in 3.usposabljanju. Narediti bi morali 4 Shapiro wilkove teste. 2. Predpostavka: SFERIČNOST (variance razlik so enake) - to preverimo s Mauchljev preizkus sferičnosti. H0:Variance razlike so enake (sferičnost izpolnjena), H1: Variance razlik so različne. Če niso enake, moramo popraviti stopinje prostosti.

H0: povprečjeZačetek=povp1=povp2=povp3

H1: vsaj eno povprečje je različno

Širok format spremenimo v dolg format - to naredimo s funkcijo pivot_longer.V cols damo spremenljivke.

library(rstatix)
## Warning: package 'rstatix' was built under R version 4.4.2
## 
## Attaching package: 'rstatix'
## The following object is masked from 'package:stats':
## 
##     filter
library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.4.2
## Warning: package 'forcats' was built under R version 4.4.2
## Warning: package 'lubridate' was built under R version 4.4.2
## ── Attaching core tidyverse packages ────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2
## ── Conflicts ──────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks rstatix::filter(), stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
podatki_dolg <- podatki %>%
                 pivot_longer(
                   cols = c("Zacetek", "PU1", "PU2", "PU3"),
                   names_to = "Obdobje", 
                   values_to = "Proizvodi") %>%
                 as.data.frame()

podatki_dolg$Obdobje <- factor(podatki_dolg$Obdobje,
                               levels = c("Zacetek", "PU1", "PU2", "PU3"),
                               labels = c("Zacetek", "PU1", "PU2", "PU3"))

head(podatki_dolg, 10)
##    ID Obdobje Proizvodi
## 1   1 Zacetek        85
## 2   1     PU1       133
## 3   1     PU2       136
## 4   1     PU3       134
## 5   2 Zacetek        90
## 6   2     PU1       123
## 7   2     PU2       138
## 8   2     PU3       135
## 9   3 Zacetek        96
## 10  3     PU1       120

S funkcijo ggboxplot narišemo porazdelitev izdelkov po obdobjih.

library(ggpubr)
## Warning: package 'ggpubr' was built under R version 4.4.2
ggboxplot(podatki_dolg, 
          x = "Obdobje", 
          y = "Proizvodi",
          xlab = "Obdobje",
          ylab = "Število proizvodov")

Največje razlike so med PU1 in Začetek.

library(rstatix)

podatki_dolg %>%
  group_by(Obdobje) %>%
  shapiro_test(Proizvodi)
## # A tibble: 4 × 4
##   Obdobje variable  statistic      p
##   <fct>   <chr>         <dbl>  <dbl>
## 1 Zacetek Proizvodi     0.975 0.545 
## 2 PU1     Proizvodi     0.959 0.179 
## 3 PU2     Proizvodi     0.957 0.146 
## 4 PU3     Proizvodi     0.947 0.0685

Normalnost je izpolnjena, ker je p večja od 0,05!

library(rstatix)

podatki_dolg %>%
  group_by(Obdobje) %>%
  get_summary_stats(Proizvodi, type = "mean_sd")
## # A tibble: 4 × 5
##   Obdobje variable      n  mean    sd
##   <fct>   <fct>     <dbl> <dbl> <dbl>
## 1 Zacetek Proizvodi    38  89.9  4.53
## 2 PU1     Proizvodi    38 125.   4.97
## 3 PU2     Proizvodi    38 135.   5.50
## 4 PU3     Proizvodi    38 138.   4.95

Če primerjamo povprečja so bili najbolj produktivni po PU3.

library(rstatix)

#rANOVA
ANOVA_results <- anova_test(dv = Proizvodi, #Dependent variable
                            wid = ID, #Subject identifier
                            within = Obdobje, #Within-subject factor variable
                            data = podatki_dolg)

ANOVA_results #Summary of results.
## ANOVA Table (type III tests)
## 
## $ANOVA
##    Effect DFn DFd      F        p p<.05   ges
## 1 Obdobje   3 111 760.45 8.09e-74     * 0.938
## 
## $`Mauchly's Test for Sphericity`
##    Effect     W     p p<.05
## 1 Obdobje 0.981 0.983      
## 
## $`Sphericity Corrections`
##    Effect   GGe       DF[GG]    p[GG] p[GG]<.05   HFe       DF[HF]
## 1 Obdobje 0.987 2.96, 109.54 7.01e-73         * 1.082 3.25, 120.11
##      p[HF] p[HF]<.05
## 1 8.09e-74         *

Preverjamo, če so povprečja vseh 4 spremenljivk enaka. Gledamo p-vrednost v prvi vrstici (8,09e-74). P-vrednost v drugi vrstici (Mauchly test) preverjamo sferičnost - v tem primeru je sferičnost je izpolnjena (p večji od 0.05). Če je sferičnost izpolnjena VEDNO gledaš prvo vrstico brez popravkov. Če bi bila kršena sferičnost (p je manjši od 0,05) gledamo prvo p.vrednost v zadnji vrstici (7.01e-73).

get_anova_table(ANOVA_results, correction = "auto")
## ANOVA Table (type III tests)
## 
##    Effect DFn DFd      F        p p<.05   ges
## 1 Obdobje   3 111 760.45 8.09e-74     * 0.938

ges pomeni eta squared - velikost učinka (0,938). Gledati moram lestvico za ES (eta kvadrat), ki se uporablja za analizo variance. 0,938 pade med velike razlike.

library(effectsize)
## 
## Attaching package: 'effectsize'
## The following objects are masked from 'package:rstatix':
## 
##     cohens_d, eta_squared
interpret_eta_squared(0.938, rules = "cohen1992")
## [1] "large"
## (Rules: cohen1992)

S pairwise_t_test ugotavljamo kje so razlike. Gledamo vse možne kombinacije. 6 kombinacij. Pri p.adj moramo pomnožiti p-vrednost s 6.

library(rstatix)

#Comparing all paires of variables
pwc <- podatki_dolg %>%
  pairwise_t_test(Proizvodi ~ Obdobje, 
                  paired = TRUE,
                  p.adjust.method = "bonferroni")

pwc
## # A tibble: 6 × 10
##   .y.      group1 group2    n1    n2 statistic    df        p    p.adj
## * <chr>    <chr>  <chr>  <int> <int>     <dbl> <dbl>    <dbl>    <dbl>
## 1 Proizvo… Zacet… PU1       38    38    -32.4     37 9.74e-29 5.84e-28
## 2 Proizvo… Zacet… PU2       38    38    -39.7     37 6.28e-32 3.77e-31
## 3 Proizvo… Zacet… PU3       38    38    -44.7     37 8.16e-34 4.9 e-33
## 4 Proizvo… PU1    PU2       38    38     -8.35    37 4.88e-10 2.93e- 9
## 5 Proizvo… PU1    PU3       38    38    -11.6     37 6.70e-14 4.02e-13
## 6 Proizvo… PU2    PU3       38    38     -2.49    37 1.7 e- 2 1.04e- 1
## # ℹ 1 more variable: p.adj.signif <chr>

Friedman ANOVA uporabimo, če normalnost ni izpolnjena. H0: lokacije porazdelitve števila proiz. izdelkov so enake. V tem primeru bi zavrnili H0 (p je manjši kot 0,05)

library(rstatix)

#Friedman ANOVA.
FriedmanANOVA <- friedman_test(Proizvodi ~ Obdobje | ID,
                               data = podatki_dolg)

FriedmanANOVA #Summary of results.
## # A tibble: 1 × 6
##   .y.           n statistic    df        p method       
## * <chr>     <int>     <dbl> <dbl>    <dbl> <chr>        
## 1 Proizvodi    38      98.9     3 2.69e-21 Friedman test

Velikost učinka pri Friedman ANOVA merimo s Kenddals_w statistiko. dobimo 0,87 kar pomeni almost perfect agreement (to pomeni, da skoraj vsi so najbolj produktivni po pu3, nato pu2, pu1, najmanj v začetek). zelo so skladni kdaj so najbolj produktivni in kdaj najmanj.

library(effectsize)
effectsize::kendalls_w(Proizvodi ~ Obdobje | ID,
                       data = podatki_dolg)
## Warning: 2 block(s) contain ties.
## Kendall's W |       95% CI
## --------------------------
## 0.87        | [0.83, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].
interpret_kendalls_w(0.87)
## [1] "almost perfect agreement"
## (Rules: landis1977)

NEPARAMETRIČNI PREIZKUS: WILCOXONOV TEST PREZNAČNIH RANGOV.

library(rstatix)

#Wilcoxon signed rank tests - comparing all possible pairs.
pairs_nonpar <- wilcox_test(Proizvodi ~ Obdobje, 
                            paired = TRUE, 
                            p.adjust.method = "bonferroni",
                            data = podatki_dolg)

pairs_nonpar
## # A tibble: 6 × 9
##   .y.       group1  group2    n1    n2 statistic            p    p.adj
## * <chr>     <chr>   <chr>  <int> <int>     <dbl>        <dbl>    <dbl>
## 1 Proizvodi Zacetek PU1       38    38        0  0.0000000794  4.76e-7
## 2 Proizvodi Zacetek PU2       38    38        0  0.000000079   4.74e-7
## 3 Proizvodi Zacetek PU3       38    38        0  0.0000000796  4.78e-7
## 4 Proizvodi PU1     PU2       38    38       24  0.000000514   3.08e-6
## 5 Proizvodi PU1     PU3       38    38        0  0.0000000794  4.76e-7
## 6 Proizvodi PU2     PU3       38    38      202. 0.041         2.45e-1
## # ℹ 1 more variable: p.adj.signif <chr>