Publikacje naukowe · Otwarta nauka

Analiza wieloświatowa (multiverse) i krzywa specyfikacji – pełny workflow w R z wnioskowaniem łącznym

Jeden zbiór danych, kilkadziesiąt rozsądnych ścieżek analizy. Pokazujemy, jak wybrać decyzje do multiverse, policzyć wszystkie specyfikacje pętlą w R, narysować dwupanelową krzywą specyfikacji i sprawdzić testem permutacyjnym, czy krzywa jako całość świadczy o efekcie.

Strona główna / Blog / Multiverse i krzywa specyfikacji
W skrócie: Analiza wieloświatowa (multiverse) powtarza tę samą analizę dla wszystkich kombinacji rozsądnych decyzji o danych i modelu, a krzywa specyfikacji pokazuje wyniki od najsilniejszego do najsłabszego razem z decyzjami, które za nimi stoją. Do multiverse wpuszczaj tylko opcje naprawdę równoważne i zawsze dodaj wnioskowanie łączne: permutację albo bootstrap przy H0, bo liczba istotnych specyfikacji sama niczego nie dowodzi.

Stopnie swobody badacza i ogród rozwidlających się ścieżek

Między surowymi danymi a jedną liczbą w tabeli wyników leży kilkanaście decyzji. Kogo wykluczyć? Co zrobić z odstającymi czasami reakcji? Jak policzyć wynik? Które kowarianty dodać? Wicherts i in. (2016) zebrali 34 takie stopnie swobody badacza, od sformułowania hipotezy po raportowanie. Simmons i in. (2011) pokazali, że nieujawniona elastyczność w tych decyzjach radykalnie zwiększa odsetek wyników fałszywie pozytywnych. Jak to wygląda w symulacji, pokazujemy w haśle p-hacking.

Gelman i Loken (2014) zwrócili uwagę na subtelniejszą wersję problemu, nazywaną ogrodem rozwidlających się ścieżek (garden of forking paths). Badacz nie musi niczego „łowić”. Wystarczy, że wykona jedną analizę, ale jej postać wybierze po zobaczeniu danych. Przy innych danych wybrałby inną ścieżkę, więc nominalne α = .05 przestaje opisywać rzeczywisty błąd I rodzaju. Że nie chodzi o teorię, pokazał projekt Silberzahna i in. (2018): 29 zespołów analizowało ten sam zbiór danych i 69% z nich uzyskało istotny efekt, a 31% nie.

Są dwie odpowiedzi na ten problem. Pierwsza to zamknięcie ścieżek przed zebraniem danych, czyli prerejestracja. Druga to pokazanie wszystkich rozsądnych ścieżek naraz. Ten wpis dotyczy drugiej: analizy wieloświatowej i krzywej specyfikacji z wnioskowaniem łącznym. Krótkie wprowadzenie do krzywej jako jednego z elementów analizy wrażliwości znajdziesz we wpisie analiza wrażliwości w artykule, a definicję w haśle krzywa specyfikacji.

Multiverse danych, multiverse modeli i krzywa specyfikacji

Steegen i in. (2016) zaproponowali analizę wieloświatową (multiverse analysis). Przetwarzanie surowych danych wymaga wyborów między kilkoma rozsądnymi sposobami wykluczania, przekształcania i kodowania. Zamiast jednego zbioru danych tworzysz więc wszystkie zbiory wynikające z tych wyborów i na każdym powtarzasz analizę. Wynik pokazuje, jak bardzo wniosek zależy od arbitralnych decyzji i które z nich ważą najwięcej. Steegen i in. skupili się na tym multiverse danych. Ten sam pomysł stosuje się do decyzji o modelu: zestawu kowariantów, rozkładu reszt czy estymatora. Mówi się wtedy o multiverse modeli. W praktyce obie przestrzenie krzyżuje się ze sobą.

Simonsohn i in. (2020) zauważają, że decyzje analityczne bywają jednocześnie uzasadnione, arbitralne i motywowane. Wprowadzają przy tym zmienność, której nie widać w błędach standardowych. Ich analiza krzywej specyfikacji (specification curve analysis) ma trzy kroki:

  1. wyznacz zbiór specyfikacji uzasadnionych teoretycznie, poprawnych statystycznie i nieredundantnych;
  2. pokaż wyniki na wykresie, z którego widać, które decyzje zmieniają wynik;
  3. przeprowadź wnioskowanie łączne ze wszystkich specyfikacji.

Trzeci krok odróżnia krzywą specyfikacji od samego wykresu wielu wariantów. Pytanie nie brzmi „ile wariantów jest istotnych?”, tylko „czy taka krzywa byłaby prawdopodobna, gdyby efektu nie było?”.

Które decyzje wpuścić do multiverse

Multiverse jest tak dobry jak lista decyzji. Del Giudice i Gangestad (2021) ostrzegają, że gdy alternatywy nie są naprawdę arbitralne, analiza może wprowadzać w błąd, a prawdziwy efekt może zginąć wśród słabo uzasadnionych wariantów. Proponują sprawdzić każdą decyzję pod kątem trzech rodzajów nierównoważności: pomiaru (czy opcje mierzą to samo), efektu (czy dotyczą tego samego efektu) oraz mocy i precyzji. Z tej oceny wynikają trzy typy decyzji:

  • typ E (równoważność): opcje są naprawdę wymienne, więc wchodzą do multiverse;
  • typ N (nierównoważność): jedna opcja jest lepsza z góry, więc nie ma jej w krzywej, tylko wybierasz ją w planie;
  • typ U (niepewność): nie wiadomo, a wtedy analizę wieloświatową prowadzisz jawnie jako eksploracyjną.

W praktyce oznacza to zasadę „arbitralne, ale rozsądne”. Każdą opcję musisz umieć obronić przed recenzentem jako tak samo dobrą jak pozostałe. Model z pominiętą znaną zmienną zakłócającą, reguła wykluczeń, której nikt w danej dziedzinie nie stosuje, albo zmiana konstruktu w zmiennej zależnej nie są alternatywami, tylko błędami. Ich dodanie sztucznie poszerza krzywą.

Przykład: krótki trening uważności a efekt Stroopa

Studenci (N = 160) zostali losowo przydzieleni do 8-minutowego ćwiczenia uważności oddechu albo do słuchania audiobooka (zob. randomizacja). Przed manipulacją i po niej wykonali zadanie Stroopa. Wynikiem jest efekt interferencji, czyli o ile milisekund wolniej odpowiadano w próbach niezgodnych niż w zgodnych. Hipoteza była kierunkowa: po ćwiczeniu uważności interferencja będzie mniejsza. Dane są symulowane (set.seed(2190)): 15 360 prób, po 48 zgodnych i 48 niezgodnych na osobę.

W danych są typowe problemy. Odpowiedzi szybsze niż 200 ms stanowią 1.33% prób, a wolniejsze niż 2000 ms 1.43%. Dziewięć osób miało trafność poniżej 80%, z czego trzy poniżej 70%. Badacz ma więc cztery decyzje do podjęcia. Dla każdej opisujemy opcje i uzasadnienie:

Cztery decyzje analityczne w przykładzie z zadaniem Stroopa, ich opcje, uzasadnienie i odrzucone warianty
DecyzjaOpcje w multiverseDlaczego są rozsądneCzego nie włączono
Odstające czasy reakcji (multiverse danych)granice 200–2000 ms; ±2.5 SD w osobie i warunku; 3 MAD w osobie i warunkustałe i względne punkty odcięcia to standardowe sposoby ograniczania wpływu odstających czasów, których skutki porównuje Ratcliff (1993); MAD jest odporną miarą rozproszenia (Leys i in., 2013)brak przycinania: antycypacje poniżej 200 ms nie są reakcją na bodziec
Wynik (multiverse danych)średnia albo mediana czasu w warunkumediana jest mniej wrażliwa na obserwacje odstające (Ratcliff, 1993)wskaźnik łączący czas i trafność, bo mierzy inny konstrukt (nierównoważność pomiaru)
Wykluczenie osób (multiverse danych)trafność ≥ 80% albo ≥ 70%oba progi odrzucają osoby odpowiadające przypadkowo (zob. nieuważne odpowiedzi)wykluczanie według wyniku interferencji, bo zależy od zmiennej zależnej
Kowarianty (multiverse modeli)brak; interferencja wyjściowa; interferencja wyjściowa i wiekprzy losowym przydziale każdy z tych modeli szacuje ten sam efekt bez obciążeniakowarianty zmierzone po manipulacji

Daje to 3 × 2 × 2 × 3 = 36 specyfikacji. Zwróć uwagę na ostatnią decyzję. Model z pomiarem wyjściowym (ANCOVA) jest przy korelacji przed–po r = .58 wyraźnie precyzyjniejszy niż porównanie samych wyników końcowych. To nierównoważność mocy w rozumieniu Del Giudice i Gangestada (2021), więc w planie analizy ANCOVA powinna być wybrana z góry (zob. ANCOVA w badaniu przed i po). Model bez kowariantów zostawiamy w krzywej celowo, żeby pokazać, jak precyzja wpływa na odsetek istotnych wyników.

Analiza zapisana w prerejestracji brzmiała: granice 200–2000 ms, średnia, trafność ≥ 80% i ANCOVA z pomiarem wyjściowym. Dała ona różnicę b = −17.6 ms, 95% CI [−29.5, −5.7], p = .004 (n = 151). Pytanie brzmi: na ile ten wynik zależy od przyjętych decyzji?

Krzywa specyfikacji krok po kroku

Obliczenia wymagają tylko bazowego R. Pakiety specr i multiverse automatyzują te kroki, ale pętla daje pełną kontrolę nad tym, co liczysz:

  1. Siatka decyzji. expand.grid() tworzy ramkę danych, w której każdy wiersz to jedna specyfikacja.
  2. Wynik osoby raz na wersję. Interferencja zależy tylko od przycinania i miary, więc liczysz ją sześć razy, a nie 36.
  3. Pętla po wierszach. Dla każdej specyfikacji wykluczasz osoby, dopasowujesz lm() i zapisujesz b, SE, 95% CI i p.
  4. Sortowanie od najsilniejszego efektu i rysunek.

Wszystkie 36 oszacowań jest ujemnych, od −19.9 do −6.5 ms, z medianą −14.5 ms. Wynik istotny (p < .05) dało 19 z 36 specyfikacji, czyli nieco ponad połowa. Rysunek 1 pokazuje, skąd ten podział.

Krzywa specyfikacji 36 analiz: wszystkie różnice ujemne od −19.9 do −6.5 ms, mediana −14.5 ms, 19 istotnych; w dolnym panelu wszystkie modele bez kowariantów leżą w nieistotnej, prawej części krzywej
Rysunek 1. Krzywa specyfikacji. Górny panel: różnica w interferencji (uważność − kontrola) z 95% CI w każdej z 36 specyfikacji; czerwone – p < .05, niebieska linia – mediana. Dolny panel: opcje, z których złożona jest dana specyfikacja. Obliczenia: R (dane symulowane).

Dolny panel czyta się kolumnami. Najwięcej zmienia decyzja o kowariantach. Bez pomiaru wyjściowego żadna z 12 specyfikacji nie jest istotna (mediana −10.2 ms). Z samym pomiarem wyjściowym istotnych jest 10 z 12 (mediana −16.7 ms), a z pomiarem wyjściowym i wiekiem 9 z 12 (−16.1 ms). Drugą decyzją jest miara: średnia dała 12 istotnych wyników z 18, mediana 7. Łagodniejszy próg trafności (≥ 70%) dał 12 istotnych z 18, a ostrzejszy 7. Reguła przycinania prawie nic nie zmienia (6–7 istotnych z 12 przy każdej regule). Wniosek dla czytelnika: kierunek efektu jest stały, a o istotności decyduje głównie precyzja modelu.

Wnioskowanie łączne: permutacja przy H0

Czy 19 istotnych wyników na 36 to dużo? Tego nie widać z samego wykresu. Specyfikacje liczone na tych samych danych są silnie skorelowane, więc nie można ich traktować jak 36 niezależnych testów. Nie pomoże też zwykła poprawka na wielokrotne porównania. Potrzebny jest rozkład całej krzywej przy prawdziwej hipotezie zerowej.

W eksperymencie z losowym przydziałem najprościej uzyskać go testem permutacyjnym. Etykiety grup tasujesz losowo, co niszczy ewentualny efekt, ale zachowuje wszystko inne: rozkłady czasów, wykluczenia, korelację z pomiarem wyjściowym. Na każdej permutacji liczysz wszystkie 36 specyfikacji i zapisujesz statystyki opisujące krzywą. W przykładzie są to mediana oszacowań i liczba istotnych specyfikacji w przewidywanym kierunku. Wartość p to odsetek permutacji, w których krzywa była co najmniej tak skrajna jak w danych.

Rozkład mediany oszacowań i liczby istotnych ujemnych specyfikacji w 1000 permutacjach przy H0; mediana z danych −14.5 ms, p = .013; 19 istotnych specyfikacji, p = .012
Rysunek 2. Wnioskowanie łączne. Szare słupki: wartości statystyk krzywej w 1000 permutacjach etykiet grup (prawdziwa H0). Czerwona linia: wartość z danych. Prawy panel w skali pierwiastkowej, bo w większości permutacji nie było żadnego istotnego wyniku. Obliczenia: R (dane symulowane).

Wyniki (testy jednostronne, zgodne z kierunkową hipotezą):

  • mediana oszacowań: −14.5 ms w danych, a 95% permutacji mieściło się w przedziale od −12.4 do 14.0 ms; p = .013;
  • liczba istotnych ujemnych specyfikacji: 19 w danych, przy H0 mediana 0, a 95. centyl 5; p = .012.

Krzywa jako całość jest więc mało prawdopodobna, gdy efektu nie ma, choć prawie połowa pojedynczych specyfikacji jest nieistotna. Permutacje pokazują też drugą stronę medalu. W 10.4% z nich choć jedna specyfikacja dała istotny wynik ujemny, choć efektu nie było. Tyle wynosiłoby ryzyko błędu badacza, który przegląda 36 wariantów i raportuje ten „udany”. Jedna z góry wybrana analiza z testem dwustronnym daje taki wynik w 2.5% przypadków.

Dane obserwacyjne. Gdy nie ma losowego przydziału, tasowanie predyktora niszczy też jego związki z kowariantami. Rozkład przy H0 uzyskuje się wtedy bootstrapem z wymuszoną hipotezą zerową. W każdej specyfikacji od zmiennej zależnej odejmujesz oszacowany efekt (y* = y − b · x), dzięki czemu efekt w tak zmienionych danych wynosi zero. Następnie losujesz obserwacje ze zwracaniem i na każdej próbie liczysz całą krzywą. Dalej postępujesz jak przy permutacjach.

Jak opisać analizę wieloświatową

Czytelnik musi zobaczyć cztery rzeczy:

  • Listę decyzji z opcjami i uzasadnieniem, najlepiej w tabeli, oraz informację, co celowo pominięto i dlaczego.
  • Podsumowanie krzywej: liczbę specyfikacji, zakres i medianę oszacowań, liczbę wyników istotnych w każdym kierunku.
  • Decyzje, które zmieniają wynik: nazwane wprost, z liczbami dla każdej opcji.
  • Wnioskowanie łączne: metodę (permutacja albo bootstrap, liczba powtórzeń), statystyki i ich wartości p.

Rysunek krzywej i pełną tabelę 36 wyników umieść w materiałach uzupełniających, razem z kodem (zob. pakiet replikacyjny; o układzie rysunków piszemy we wpisie tabele i rysunki do artykułu). W tekście wystarczy akapit:

Aby ocenić, na ile wynik zależy od decyzji analitycznych, przeprowadzono analizę krzywej specyfikacji (Simonsohn i in., 2020). Skrzyżowano cztery decyzje: regułę usuwania odstających czasów reakcji (granice 200–2000 ms, ±2.5 SD lub 3 MAD w osobie i warunku), miarę interferencji (średnia lub mediana), próg wykluczenia osób (trafność co najmniej 80% lub 70%) oraz zestaw kowariantów (brak; interferencja wyjściowa; interferencja wyjściowa i wiek), co dało 36 specyfikacji. Wszystkie oszacowania różnicy między grupą uważności a grupą kontrolną były ujemne (od −19.9 do −6.5 ms, mediana −14.5 ms), a 19 z 36 osiągnęło p < .05. Żadna z 12 specyfikacji bez pomiaru wyjściowego nie była istotna, natomiast istotnych było 19 z 24 specyfikacji z tym kowariantem. W teście permutacyjnym (1000 permutacji etykiet grup) mediana oszacowań, p = .013, i liczba istotnych ujemnych specyfikacji, p = .012, były większe niż przy prawdziwej hipotezie zerowej (testy jednostronne). Analiza zaplanowana w prerejestracji dała b = −17.6 ms, 95% CI [−29.5, −5.7], p = .004. Pełną listę specyfikacji przedstawiono w tabeli S1 i na rysunku S1.

Gdy recenzent prosi o taką analizę po złożeniu artykułu, przydadzą się wzory z wpisu odpowiedź recenzentom. Jeśli chcesz, żeby ktoś z zewnątrz sprawdził listę decyzji i kod permutacji, zrobimy to w ramach wsparcia statystycznego dla badaczy.

Multiverse a prerejestracja

Prerejestracja i multiverse nie konkurują ze sobą. Nosek i in. (2018) wskazują, że główną funkcją prerejestracji jest odróżnienie testowania przewidywań od wyjaśniania wyników po fakcie. Analiza główna zapisana z góry zachowuje swoje znaczenie konfirmacyjne. Krzywa specyfikacji odpowiada na inne pytanie: czy wniosek przetrwa decyzje, które równie dobrze mogły zapaść inaczej. Najlepiej zarejestrować obie rzeczy: analizę główną oraz listę decyzji do multiverse z wybranymi statystykami łącznymi. Jak przygotować taki plan, opisujemy we wpisie plan analizy statystycznej i prerejestracja.

Bez prerejestracji multiverse jest szczególnie cenny, bo pokazuje ścieżki, których autor nie wybrał. Pozostaje jednak analizą, którą trzeba opisać uczciwie. Jeśli lista decyzji powstała po obejrzeniu danych, napisz to. Gdy decyzje są typu U, nazwij całą analizę eksploracyjną (Del Giudice i Gangestad, 2021). Przy reanalizie cudzego artykułu multiverse służy do oceny kruchości opublikowanego wyniku, co pokazujemy we wpisie reanaliza i wtórna analiza danych. Krzywa specyfikacji pasuje też do formatu registered report i do praktyk otwartej nauki, bo w obu przypadkach liczy się pełny obraz, a nie jedna ścieżka.

Przykład reprodukowalny

Przykład w R

Cztery przepisy uruchamiane po kolei w jednej sesji. Obliczenia (przepisy 1–3) wymagają tylko bazowego R: expand.grid(), pętli i lm(). Permutacje to 36 000 modeli i trwają kilkanaście sekund. Do rysunku potrzebne są pakiety ggplot2 i patchwork.

1. Dane i wynik osoby

Dane wejściowe: 160 osób, losowy przydział do ćwiczenia uważności albo kontroli, po 48 prób zgodnych i 48 niezgodnych w zadaniu Stroopa, interferencja przed manipulacją i wiek (dane symulowane, set.seed(2190)).

Co odczytać: liczbę osób i prób, liczbę osób z trafnością poniżej 80% i 70%, odsetek prób poniżej 200 ms i powyżej 2000 ms oraz odsetek poprawnych prób usuniętych przez każdą z trzech reguł przycinania. Funkcja interferencja() liczy wynik jednej osoby dla wybranej reguły i miary.

R
# Dane symulowane: eksperyment z losowym przydziałem. 8-minutowe ćwiczenie uważności oddechu albo słuchanie audiobooka
# (kontrola), potem zadanie Stroopa: 48 prób zgodnych i 48 niezgodnych na osobę. Wynik: efekt interferencji (ms).
set.seed(2190)
n_os <- 160
osoby <- data.frame(id = 1:n_os, grupa = factor(rep(c("kontrola", "uwaznosc"), each = n_os / 2),
                                                 levels = c("kontrola", "uwaznosc")),
                    wiek = round(runif(n_os, 19, 35)))
tau <- rnorm(n_os, 85, 30)                                    # prawdziwa interferencja osoby (ms)
osoby$interf_pre <- round(tau + rnorm(n_os, 0, 25))           # interferencja zmierzona przed manipulacją
tau_po <- tau - 16 * (osoby$grupa == "uwaznosc") + rnorm(n_os, 0, 12)
baza <- 640 + 4 * (osoby$wiek - 25) + rnorm(n_os, 0, 80)      # ogólne tempo reakcji osoby
nieuwazni <- c(sample(1:80, 4), sample(81:160, 5))   # osoby odpowiadające bez skupienia
trafnosc_os <- ifelse(seq_len(n_os) %in% nieuwazni, runif(n_os, .62, .78), runif(n_os, .90, .99))

proby <- expand.grid(proba = 1:48, warunek = c("zgodna", "niezgodna"), id = 1:n_os)
proby <- proby[, c("id", "warunek", "proba")]
i <- proby$id
proby$rt <- baza[i] + tau_po[i] * (proby$warunek == "niezgodna") + rnorm(nrow(proby), 0, 90) + rexp(nrow(proby), 1 / 120)
los <- runif(nrow(proby))
proby$rt[los < .015] <- runif(sum(los < .015), 100, 200)                     # antycypacje
proby$rt[los > .985] <- runif(sum(los > .985), 2000, 4000)                   # zawieszenia uwagi
proby$rt[i %in% nieuwazni] <- runif(sum(i %in% nieuwazni), 250, 900)         # szybkie, przypadkowe odpowiedzi
proby$rt <- round(proby$rt)
proby$poprawna <- runif(nrow(proby)) < trafnosc_os[i]

# Trafność osoby (odsetek poprawnych odpowiedzi) i opis danych
osoby$trafnosc <- as.numeric(tapply(proby$poprawna, proby$id, mean))
c(osob = n_os, prob = nrow(proby), trafnosc_ponizej_80 = sum(osoby$trafnosc < .80), trafnosc_ponizej_70 = sum(osoby$trafnosc < .70))
round(c(rt_ponizej_200 = mean(proby$rt < 200), rt_powyzej_2000 = mean(proby$rt > 2000)) * 100, 2)   # odsetek prób (%)

# Reguła przycinania czasów reakcji w obrębie osoby (SD i MAD liczone osobno dla każdego warunku)
zostaw_proby <- function(rt, warunek, przycinanie) {
  switch(przycinanie,
    "200–2000 ms" = rt >= 200 & rt <= 2000,
    "±2.5 SD" = abs(rt - ave(rt, warunek)) <= 2.5 * ave(rt, warunek, FUN = sd),
    "3 MAD" = abs(rt - ave(rt, warunek, FUN = median)) <= 3 * ave(rt, warunek, FUN = mad))
}
# Interferencja jednej osoby: niezgodne − zgodne, tylko poprawne próby; miara: średnia albo mediana
interferencja <- function(rt, warunek, przycinanie, miara) {
  zostaw <- zostaw_proby(rt, warunek, przycinanie)
  f <- if (miara == "średnia") mean else median
  f(rt[zostaw & warunek == "niezgodna"]) - f(rt[zostaw & warunek == "zgodna"])
}
poprawne <- proby[proby$poprawna, ]
# Odsetek poprawnych prób usuniętych przez każdą regułę (%)
sapply(c("200–2000 ms", "±2.5 SD", "3 MAD"), function(r)
  round(100 * mean(!unlist(lapply(split(poprawne, poprawne$id), function(d) zostaw_proby(d$rt, d$warunek, r)))), 1))
Wynik w R
               osob                prob trafnosc_ponizej_80 trafnosc_ponizej_70 
                160               15360                   9                   3 
 rt_ponizej_200 rt_powyzej_2000 
           1.33            1.43 
200–2000 ms     ±2.5 SD       3 MAD 
        2.8         2.8         5.1 

2. Siatka decyzji, pętla i podsumowanie krzywej

Dane wejściowe: obiekty z przepisu 1; cztery decyzje o 3, 2, 2 i 3 opcjach.

Co odczytać: liczbę specyfikacji (36), opis interferencji w grupach i jej korelację z pomiarem wyjściowym, liczbę ujemnych i istotnych oszacowań, medianę i zakres b, a w czterech małych tabelach medianę b i liczbę istotnych wyników dla każdej opcji. Na końcu pełna tabela 36 specyfikacji od najsilniejszego efektu.

R
# Krzywa specyfikacji: wszystkie kombinacje czterech decyzji (expand.grid) i jedna pętla
decyzje <- expand.grid(
  przycinanie = c("200–2000 ms", "±2.5 SD", "3 MAD"),                       # reguła dla odstających czasów reakcji
  miara = c("średnia", "mediana"),                                           # definicja wyniku
  prog_trafnosci = c(.80, .70),                                              # wykluczenie osób
  kowarianty = c("brak", "pomiar wyjściowy", "pomiar wyjściowy + wiek"),     # model
  stringsAsFactors = FALSE)
nrow(decyzje)

# Wynik osoby zależy tylko od przycinania i miary (6 wersji) – policz go raz dla każdej wersji
wersje <- unique(decyzje[c("przycinanie", "miara")])
wyniki_osob <- list()
for (w in seq_len(nrow(wersje))) {
  y <- sapply(split(poprawne, poprawne$id), function(d)
    interferencja(d$rt, d$warunek, wersje$przycinanie[w], wersje$miara[w]))
  wyniki_osob[[paste(wersje$przycinanie[w], wersje$miara[w])]] <- y
}

# Opis: interferencja (200–2000 ms, średnia) u osób z trafnością ≥ 80% i jej korelacja z pomiarem wyjściowym
y1 <- wyniki_osob[["200–2000 ms średnia"]][osoby$trafnosc >= .80]
g1 <- osoby$grupa[osoby$trafnosc >= .80]
table(g1)
round(rbind(M = tapply(y1, g1, mean), SD = tapply(y1, g1, sd)), 1)
round(cor(y1, osoby$interf_pre[osoby$trafnosc >= .80]), 2)

wzory <- c("brak" = "y ~ grupa", "pomiar wyjściowy" = "y ~ grupa + interf_pre",
           "pomiar wyjściowy + wiek" = "y ~ grupa + interf_pre + wiek")

# Jedna specyfikacja: różnica uważność − kontrola w interferencji (ms), SE, 95% CI, p
licz_spec <- function(s, grupa) {
  d <- osoby
  d$grupa <- grupa
  d$y <- wyniki_osob[[paste(decyzje$przycinanie[s], decyzje$miara[s])]]
  d <- d[d$trafnosc >= decyzje$prog_trafnosci[s], ]
  fit <- lm(as.formula(wzory[[decyzje$kowarianty[s]]]), data = d)
  m <- summary(fit)$coefficients["grupauwaznosc", ]
  c(n = nrow(d), b = m[[1]], se = m[[2]], t = m[[3]], p = m[[4]],
    dolny = m[[1]] - qt(.975, fit$df.residual) * m[[2]], gorny = m[[1]] + qt(.975, fit$df.residual) * m[[2]])
}

wyn <- matrix(NA, nrow(decyzje), 7, dimnames = list(NULL, c("n", "b", "se", "t", "p", "dolny", "gorny")))
for (s in seq_len(nrow(decyzje))) wyn[s, ] <- licz_spec(s, osoby$grupa)
spec <- cbind(decyzje, wyn)
spec <- spec[order(spec$b), ]
spec$ranga <- seq_len(nrow(spec))

# Podsumowanie krzywej
c(specyfikacji = nrow(spec), ujemnych = sum(spec$b < 0), istotnych_ujemnych = sum(spec$b < 0 & spec$p < .05),
  istotnych_dodatnich = sum(spec$b > 0 & spec$p < .05))
round(c(mediana_b = median(spec$b), min_b = min(spec$b), max_b = max(spec$b)), 1)
range(spec$n)

# Które decyzje ważą najwięcej: mediana b i liczba istotnych specyfikacji w obrębie każdej opcji
for (k in c("przycinanie", "miara", "prog_trafnosci", "kowarianty")) {
  print(do.call(rbind, lapply(split(spec, spec[[k]]), function(d)
    data.frame(specyfikacji = nrow(d), mediana_b = round(median(d$b), 1), istotnych = sum(d$p < .05)))))
}

# Pełna tabela (do materiałów uzupełniających)
print(transform(spec[c("ranga", "przycinanie", "miara", "prog_trafnosci", "kowarianty", "n")],
                b = round(spec$b, 1), dolny = round(spec$dolny, 1), gorny = round(spec$gorny, 1), p = round(spec$p, 3)),
      row.names = FALSE)
Wynik w R
[1] 36
g1
kontrola uwaznosc 
      76       75 
   kontrola uwaznosc
M      87.2     75.7
SD     43.0     49.1
[1] 0.58
       specyfikacji            ujemnych  istotnych_ujemnych istotnych_dodatnich 
                 36                  36                  19                   0 
mediana_b     min_b     max_b 
    -14.5     -19.9      -6.5 
[1] 151 157
            specyfikacji mediana_b istotnych
±2.5 SD               12     -14.8         6
200–2000 ms           12     -14.9         7
3 MAD                 12     -14.1         6
        specyfikacji mediana_b istotnych
mediana           18     -13.2         7
średnia           18     -16.8        12
    specyfikacji mediana_b istotnych
0.7           18     -16.1        12
0.8           18     -13.2         7
                        specyfikacji mediana_b istotnych
brak                              12     -10.2         0
pomiar wyjściowy                  12     -16.7        10
pomiar wyjściowy + wiek           12     -16.1         9
 ranga przycinanie   miara prog_trafnosci              kowarianty   n     b
     1 200–2000 ms średnia            0.7        pomiar wyjściowy 157 -19.9
     2 200–2000 ms średnia            0.7 pomiar wyjściowy + wiek 157 -19.6
     3     ±2.5 SD średnia            0.7        pomiar wyjściowy 157 -19.2
     4     ±2.5 SD średnia            0.7 pomiar wyjściowy + wiek 157 -18.9
     5       3 MAD średnia            0.7        pomiar wyjściowy 157 -18.1
     6       3 MAD średnia            0.7 pomiar wyjściowy + wiek 157 -17.8
     7 200–2000 ms średnia            0.8        pomiar wyjściowy 151 -17.6
     8 200–2000 ms średnia            0.8 pomiar wyjściowy + wiek 151 -17.5
     9     ±2.5 SD średnia            0.8        pomiar wyjściowy 151 -16.9
    10 200–2000 ms mediana            0.7        pomiar wyjściowy 157 -16.7
    11     ±2.5 SD średnia            0.8 pomiar wyjściowy + wiek 151 -16.7
    12     ±2.5 SD mediana            0.7        pomiar wyjściowy 157 -16.6
    13 200–2000 ms mediana            0.7 pomiar wyjściowy + wiek 157 -16.2
    14     ±2.5 SD mediana            0.7 pomiar wyjściowy + wiek 157 -16.0
    15       3 MAD mediana            0.7        pomiar wyjściowy 157 -16.0
    16       3 MAD średnia            0.8        pomiar wyjściowy 151 -15.8
    17       3 MAD średnia            0.8 pomiar wyjściowy + wiek 151 -15.6
    18       3 MAD mediana            0.7 pomiar wyjściowy + wiek 157 -15.4
    19 200–2000 ms mediana            0.8        pomiar wyjściowy 151 -13.6
    20 200–2000 ms średnia            0.7                    brak 157 -13.6
    21     ±2.5 SD mediana            0.8        pomiar wyjściowy 151 -13.5
    22 200–2000 ms mediana            0.8 pomiar wyjściowy + wiek 151 -13.3
    23     ±2.5 SD mediana            0.8 pomiar wyjściowy + wiek 151 -13.2
    24       3 MAD mediana            0.8        pomiar wyjściowy 151 -12.8
    25     ±2.5 SD średnia            0.7                    brak 157 -12.5
    26       3 MAD mediana            0.8 pomiar wyjściowy + wiek 151 -12.5
    27       3 MAD średnia            0.7                    brak 157 -11.7
    28 200–2000 ms średnia            0.8                    brak 151 -11.5
    29     ±2.5 SD średnia            0.8                    brak 151 -10.4
    30 200–2000 ms mediana            0.7                    brak 157 -10.2
    31     ±2.5 SD mediana            0.7                    brak 157 -10.2
    32       3 MAD średnia            0.8                    brak 151  -9.6
    33       3 MAD mediana            0.7                    brak 157  -9.4
    34     ±2.5 SD mediana            0.8                    brak 151  -7.2
    35 200–2000 ms mediana            0.8                    brak 151  -7.2
    36       3 MAD mediana            0.8                    brak 151  -6.5
 dolny gorny     p
 -32.4  -7.4 0.002
 -32.1  -7.1 0.002
 -33.3  -5.1 0.008
 -33.0  -4.8 0.009
 -31.3  -5.0 0.007
 -30.9  -4.7 0.008
 -29.5  -5.7 0.004
 -29.4  -5.5 0.004
 -30.7  -3.2 0.016
 -31.6  -1.9 0.027
 -30.5  -2.9 0.018
 -31.9  -1.3 0.033
 -30.9  -1.4 0.032
 -31.2  -0.8 0.039
 -31.4  -0.5 0.043
 -28.5  -3.2 0.015
 -28.3  -2.9 0.016
 -30.7   0.0 0.050
 -27.1  -0.1 0.048
 -28.9   1.8 0.083
 -27.6   0.6 0.060
 -26.8   0.2 0.054
 -27.3   0.9 0.067
 -27.0   1.4 0.078
 -29.5   4.4 0.145
 -26.8   1.8 0.086
 -27.6   4.2 0.148
 -26.3   3.4 0.129
 -27.1   6.2 0.216
 -27.6   7.2 0.250
 -27.9   7.6 0.260
 -25.1   5.9 0.223
 -27.3   8.5 0.300
 -23.9   9.5 0.394
 -23.5   9.1 0.383
 -23.3  10.4 0.450

3. Wnioskowanie łączne: permutacje

Dane wejściowe: obiekty z przepisów 1 i 2; 1000 permutacji etykiet grup (set.seed(2191)).

Co odczytać: statystyki krzywej w danych (mediana b, liczba istotnych ujemnych specyfikacji), jednostronne p dla obu statystyk, 50. i 95. centyl liczby istotnych specyfikacji przy H0, odsetek permutacji z co najmniej jednym istotnym wynikiem ujemnym i 95% zakres mediany przy H0.

R
# Wnioskowanie łączne: rozkład statystyk krzywej przy H0 z permutacji etykiet grup (losowy przydział w eksperymencie)
statystyki <- function(grupa) {
  w <- t(sapply(seq_len(nrow(decyzje)), licz_spec, grupa = grupa))
  c(mediana_b = median(w[, "b"]), istotnych_ujemnych = sum(w[, "b"] < 0 & w[, "p"] < .05))
}
obserwowane <- statystyki(osoby$grupa)
obserwowane

set.seed(2191)
B <- 1000
h0 <- t(replicate(B, statystyki(sample(osoby$grupa))))   # 1000 × 36 modeli, kilkanaście sekund

# p łączne: odsetek permutacji co najmniej tak skrajnych jak dane (z poprawką +1)
p_mediana <- (sum(h0[, "mediana_b"] <= obserwowane[["mediana_b"]]) + 1) / (B + 1)
p_istotne <- (sum(h0[, "istotnych_ujemnych"] >= obserwowane[["istotnych_ujemnych"]]) + 1) / (B + 1)
round(c(p_mediana = p_mediana, p_istotne = p_istotne), 3)

# Jak wygląda krzywa, gdy efektu nie ma: typowa liczba istotnych specyfikacji przy H0
quantile(h0[, "istotnych_ujemnych"], c(.50, .95), type = 1)
round(mean(h0[, "istotnych_ujemnych"] > 0), 3)               # odsetek permutacji z choć jedną istotną ujemną
round(quantile(h0[, "mediana_b"], c(.025, .975)), 1)
Wynik w R
         mediana_b istotnych_ujemnych 
         -14.48854           19.00000 
p_mediana p_istotne 
    0.013     0.012 
50% 95% 
  0   5 
[1] 0.104
 2.5% 97.5% 
-12.4  14.0 

4. Dwupanelowy rysunek w ggplot2

Dane wejściowe: obiekt spec z przepisu 2.

Co odczytać: wymiary danych długich do dolnego panelu (36 specyfikacji × 4 decyzje = 144 wiersze); rysunek trafia do pliku krzywa-specyfikacji.png. Rysunek 1 w tekście powstał z tego samego kodu z motywem graficznym serwisu.

R
# Dwupanelowa krzywa specyfikacji w ggplot2: u góry oszacowania z 95% CI, na dole decyzje każdej specyfikacji
library(ggplot2)
library(patchwork)
spec$wynik <- ifelse(spec$p < .05, "p < .05", "p ≥ .05")
kolory <- c("p < .05" = "#e50101", "p ≥ .05" = "grey55")

gora <- ggplot(spec, aes(ranga, b, colour = wynik)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_hline(yintercept = median(spec$b), colour = "#1f78b4") +
  geom_errorbar(aes(ymin = dolny, ymax = gorny), width = 0) +
  geom_point(size = 2) +
  scale_colour_manual(values = kolory, name = NULL) +
  labs(x = NULL, y = "Uważność − kontrola (ms)") +
  theme_minimal() + theme(axis.text.x = element_blank())

# Dane długie do dolnego panelu: jedna linia = specyfikacja × decyzja
etykiety <- c(przycinanie = "odstające RT", miara = "wynik", prog_trafnosci = "trafność", kowarianty = "kowarianty")
dl <- do.call(rbind, lapply(names(etykiety), function(k)
  data.frame(ranga = spec$ranga, decyzja = etykiety[[k]], opcja = as.character(spec[[k]]), wynik = spec$wynik)))
dl$decyzja <- factor(dl$decyzja, levels = etykiety)
dim(dl)

dol <- ggplot(dl, aes(ranga, opcja, colour = wynik)) +
  geom_point(shape = 15) +
  facet_grid(decyzja ~ ., scales = "free_y", space = "free_y") +
  scale_colour_manual(values = kolory, guide = "none") +
  labs(x = "Specyfikacje od najsilniejszego efektu", y = NULL) +
  theme_minimal()

wykres <- gora / dol + plot_layout(heights = c(1.2, 1))
ggsave("krzywa-specyfikacji.png", wykres, width = 10, height = 7.5, dpi = 150, bg = "white")
Wynik w R
[1] 144   4

Bibliografia oprogramowania

Wyniki na tej stronie pochodzą z uruchomienia pokazanego kodu w R 4.5.3 (funkcje expand.grid(), lm(), replicate(), sample()); rysunki przygotowano w pakietach ggplot2 4.0.2 i patchwork 1.3.2.

  1. R: R Core Team. (2026). R: A language and environment for statistical computing (Version 4.5.3) [Computer software]. R Foundation for Statistical Computing. https://www.R-project.org/
  2. ggplot2 (wersja 4.0.2): Wickham, H. (2016). ggplot2: Elegant graphics for data analysis. Springer. https://ggplot2.tidyverse.org
  3. patchwork (wersja 1.3.2): Pedersen, T. L. (2025). patchwork: The composer of plots [R package]. https://doi.org/10.32614/CRAN.package.patchwork

Najczęstsze błędy

  • Opcje, których nikt by nie wybrał. Specyfikacje bez uzasadnienia rozmywają efekt albo go pompują. Każdą opcję uzasadnij tak, jakby miała być analizą główną.
  • Zmiana konstruktu nazwana „definicją wyniku”. Inna zmienna zależna to inne pytanie, a nie kolejna gałąź multiverse.
  • Liczenie istotnych specyfikacji bez wnioskowania łącznego. „19 z 36” nic nie mówi bez rozkładu przy H0. Specyfikacje są zależne, więc nie wystarczy porównać tej liczby z 5%.
  • Wykluczenia zależne od zmiennej zależnej (np. usuwanie osób o skrajnej interferencji). Takie reguły łamią logikę permutacji i otwierają nowe ścieżki.
  • Wybór „najlepszej” specyfikacji z krzywej jako wyniku głównego. To p-hacking z wykresem w tle; w przykładzie sięgnięcie po najkorzystniejszą z 36 ścieżek dawało fałszywy wynik w 10.4% permutacji.
  • Brak nazwanej decyzji kluczowej. Ogólnik „wyniki były odporne” zastąp zdaniem, od której decyzji zależy istotność, z liczbami dla każdej opcji.

Najczęstsze pytania

Czym różni się analiza multiverse od krzywej specyfikacji?

Obie liczą tę samą analizę dla wielu kombinacji decyzji. Analiza wieloświatowa (Steegen i in., 2016) kładzie nacisk na zbiory danych powstające z różnych decyzji o wykluczeniach, przekształceniach i kodowaniu. Krzywa specyfikacji (Simonsohn i in., 2020) obejmuje wszystkie decyzje analityczne, także o modelu, pokazuje na wykresie, które z nich zmieniają wynik, i dodaje wnioskowanie łączne ze wszystkich specyfikacji.

Ile specyfikacji powinna mieć analiza wieloświatowa?

Tyle, ile wynika z rozsądnych decyzji, a nie z chęci uzyskania dużej liczby. Kilka decyzji o 2–3 opcjach daje zwykle kilkadziesiąt specyfikacji, jak 36 w przykładzie. Ważniejsze od liczby jest to, by każdą opcję dało się obronić jako równie dobrą. Słabo uzasadnione warianty tylko poszerzają krzywą (Del Giudice i Gangestad, 2021).

Jak policzyć łączną wartość p dla krzywej specyfikacji?

Wybierz statystyki opisujące całą krzywą, np. medianę oszacowań i liczbę istotnych wyników w przewidywanym kierunku. Potem wygeneruj ich rozkład przy prawdziwej H0: w eksperymencie permutując etykiety grup, w danych obserwacyjnych bootstrapem po odjęciu oszacowanego efektu. Wartość p to odsetek powtórzeń z wynikiem co najmniej tak skrajnym jak w danych.

Czy analizę multiverse trzeba zaplanować w prerejestracji?

Nie jest to warunek konieczny, ale bardzo ją wzmacnia. Zarejestruj analizę główną, listę decyzji z opcjami i statystyki łączne. Jeśli lista powstała po obejrzeniu danych albo decyzje są niepewne, opisz analizę jako eksploracyjną. Prerejestracja i multiverse odpowiadają na różne pytania, więc najlepiej stosować obie.

Jak zrobić krzywą specyfikacji w R bez pakietu specr?

Utwórz siatkę decyzji funkcją expand.grid(), w pętli dopasuj model dla każdego wiersza i zapisz oszacowanie, przedział ufności i p. Posortuj wyniki i narysuj dwa panele w ggplot2: oszacowania z przedziałami oraz macierz punktów pokazującą opcje każdej specyfikacji; połącz je pakietem patchwork. Pełny kod jest w sekcji z przykładem.

Powiązane hasła

Przeczytaj też

Recenzent pyta, czy wynik przetrwa inne decyzje analityczne?

Pomożemy dobrać uzasadnione specyfikacje, policzyć krzywą z wnioskowaniem łącznym i wyjaśnić wyniki tak, żebyś mógł(a) je samodzielnie opisać i obronić.

Oferta dla badaczyPorozmawiaj o projekcie

Bibliografia i dalsza literatura naukowa (APA 7)

  1. Del Giudice, M., & Gangestad, S. W. (2021). A traveler’s guide to the multiverse: Promises, pitfalls, and a framework for the evaluation of analytic decisions. Advances in Methods and Practices in Psychological Science, 4(1), Article 2515245920954925. https://doi.org/10.1177/2515245920954925 ● Otwarty dostęp
  2. Gelman, A., & Loken, E. (2014). The statistical crisis in science. American Scientist, 102(6), 460. https://doi.org/10.1511/2014.111.460
  3. Leys, C., Ley, C., Klein, O., Bernard, P., & Licata, L. (2013). Detecting outliers: Do not use standard deviation around the mean, use absolute deviation around the median. Journal of Experimental Social Psychology, 49(4), 764–766. https://doi.org/10.1016/j.jesp.2013.03.013 ● Otwarty dostęp
  4. Nosek, B. A., Ebersole, C. R., DeHaven, A. C., & Mellor, D. T. (2018). The preregistration revolution. Proceedings of the National Academy of Sciences, 115(11), 2600–2606. https://doi.org/10.1073/pnas.1708274114 ● Otwarty dostęp: pełny tekst
  5. Ratcliff, R. (1993). Methods for dealing with reaction time outliers. Psychological Bulletin, 114(3), 510–532. https://doi.org/10.1037/0033-2909.114.3.510
  6. Silberzahn, R., Uhlmann, E. L., Martin, D. P., Anselmi, P., Aust, F., Awtrey, E., Bahník, Š., Bai, F., Bannard, C., Bonnier, E., Carlsson, R., Cheung, F., Christensen, G., Clay, R., Craig, M. A., Dalla Rosa, A., Dam, L., Evans, M. H., Flores Cervantes, I., . . . Nosek, B. A. (2018). Many Analysts, One Data Set: Making Transparent How Variations in Analytic Choices Affect Results. Advances in Methods and Practices in Psychological Science, 1(3), 337–356. https://doi.org/10.1177/2515245917747646 ● Otwarty dostęp
  7. Simmons, J. P., Nelson, L. D., & Simonsohn, U. (2011). False-positive psychology: Undisclosed flexibility in data collection and analysis allows presenting anything as significant. Psychological Science, 22(11), 1359–1366. https://doi.org/10.1177/0956797611417632
  8. Simonsohn, U., Simmons, J. P., & Nelson, L. D. (2020). Specification curve analysis. Nature Human Behaviour, 4(11), 1208–1214. https://doi.org/10.1038/s41562-020-0912-z
  9. Steegen, S., Tuerlinckx, F., Gelman, A., & Vanpaemel, W. (2016). Increasing transparency through a multiverse analysis. Perspectives on Psychological Science, 11(5), 702–712. https://doi.org/10.1177/1745691616658637
  10. Wicherts, J. M., Veldkamp, C. L. S., Augusteijn, H. E. M., Bakker, M., van Aert, R. C. M., & van Assen, M. A. L. M. (2016). Degrees of freedom in planning, running, analyzing, and reporting psychological studies: A checklist to avoid p-hacking. Frontiers in Psychology, 7, Article 1832. https://doi.org/10.3389/fpsyg.2016.01832 ● Otwarty dostęp

Cytowania BibTeX

@article{delgiudice2021traveler, author={Del Giudice, Marco and Gangestad, Steven W.}, title={{A traveler’s guide to the multiverse: Promises, pitfalls, and a framework for the evaluation of analytic decisions}}, journal={Advances in Methods and Practices in Psychological Science}, year={2021}, volume={4}, number={1}, eid={2515245920954925}, doi={10.1177/2515245920954925}}
@article{gelman2014statistical, author={Gelman, Andrew and Loken, Eric}, title={{The statistical crisis in science}}, journal={American Scientist}, year={2014}, volume={102}, number={6}, pages={460}, doi={10.1511/2014.111.460}}
@article{leys2013detecting, author={Leys, Christophe and Ley, Christophe and Klein, Olivier and Bernard, Philippe and Licata, Laurent}, title={{Detecting outliers: Do not use standard deviation around the mean, use absolute deviation around the median}}, journal={Journal of Experimental Social Psychology}, year={2013}, volume={49}, number={4}, pages={764--766}, doi={10.1016/j.jesp.2013.03.013}}
@article{nosek2018preregistration, author={Nosek, Brian A. and Ebersole, Charles R. and DeHaven, Alexander C. and Mellor, David T.}, title={{The preregistration revolution}}, journal={Proceedings of the National Academy of Sciences}, year={2018}, volume={115}, number={11}, pages={2600--2606}, doi={10.1073/pnas.1708274114}}
@article{ratcliff1993methods, author={Ratcliff, Roger}, title={{Methods for dealing with reaction time outliers}}, journal={Psychological Bulletin}, year={1993}, volume={114}, number={3}, pages={510--532}, doi={10.1037/0033-2909.114.3.510}}
@article{silberzahn2018many, author={Silberzahn, R. and Uhlmann, E. L. and Martin, D. P. and Anselmi, P. and Aust, F. and Awtrey, E. and Bahník, Š. and Bai, F. and Bannard, C. and Bonnier, E. and Carlsson, R. and Cheung, F. and Christensen, G. and Clay, R. and Craig, M. A. and Dalla Rosa, A. and Dam, L. and Evans, M. H. and Flores Cervantes, I. and Fong, N. and Gamez-Djokic, M. and Glenz, A. and Gordon-McKeon, S. and Heaton, T. J. and Hederos, K. and Heene, M. and Hofelich Mohr, A. J. and Högden, F. and Hui, K. and Johannesson, M. and Kalodimos, J. and Kaszubowski, E. and Kennedy, D. M. and Lei, R. and Lindsay, T. A. and Liverani, S. and Madan, C. R. and Molden, D. and Molleman, E. and Morey, R. D. and Mulder, L. B. and Nijstad, B. R. and Pope, N. G. and Pope, B. and Prenoveau, J. M. and Rink, F. and Robusto, E. and Roderique, H. and Sandberg, A. and Schlüter, E. and Schönbrodt, F. D. and Sherman, M. F. and Sommer, S. A. and Sotak, K. and Spain, S. and Spörlein, C. and Stafford, T. and Stefanutti, L. and Tauber, S. and Ullrich, J. and Vianello, M. and Wagenmakers, E.-J. and Witkowiak, M. and Yoon, S. and Nosek, B. A.}, title={{Many Analysts, One Data Set: Making Transparent How Variations in Analytic Choices Affect Results}}, journal={Advances in Methods and Practices in Psychological Science}, year={2018}, volume={1}, number={3}, pages={337--356}, doi={10.1177/2515245917747646}}
@article{simmons2011false, author={Simmons, Joseph P. and Nelson, Leif D. and Simonsohn, Uri}, title={{False-positive psychology: Undisclosed flexibility in data collection and analysis allows presenting anything as significant}}, journal={Psychological Science}, year={2011}, volume={22}, number={11}, pages={1359--1366}, doi={10.1177/0956797611417632}}
@article{simonsohn2020specification, author={Simonsohn, Uri and Simmons, Joseph P. and Nelson, Leif D.}, title={{Specification curve analysis}}, journal={Nature Human Behaviour}, year={2020}, volume={4}, number={11}, pages={1208--1214}, doi={10.1038/s41562-020-0912-z}}
@article{steegen2016increasing, author={Steegen, Sara and Tuerlinckx, Francis and Gelman, Andrew and Vanpaemel, Wolf}, title={{Increasing transparency through a multiverse analysis}}, journal={Perspectives on Psychological Science}, year={2016}, volume={11}, number={5}, pages={702--712}, doi={10.1177/1745691616658637}}
@article{wicherts2016degrees, author={Wicherts, Jelte M. and Veldkamp, Coosje L. S. and Augusteijn, Hilde E. M. and Bakker, Marjan and van Aert, Robbie C. M. and van Assen, Marcel A. L. M.}, title={{Degrees of freedom in planning, running, analyzing, and reporting psychological studies: A checklist to avoid p-hacking}}, journal={Frontiers in Psychology}, year={2016}, volume={7}, eid={1832}, doi={10.3389/fpsyg.2016.01832}}
R