By Leave a Comment on Modeluri de mediere

Tabelă de materii

Actualizat la 04 august 2020. Codul poate fi descărcat aici.

Introducere

În unele situații putem lua în considerare efectul indirect al unei variabile asupra unui rezultat sau a unui efect. De exemplu, condițiile de trai precare de acasă în copilărie pot diminua rezultatele învățării la școală, care ulterior au un efect negativ asupra calității vieții ulterioare, de exemplu, asupra veniturilor obținute pe parcursul vieții. Într-un alt caz, am putea lua în considerare o singură variabilă colectată în mai multe momente de timp, astfel încât să existe un efect al variabilei la momentul 1 asupra momentului 2 și al momentului 2 asupra momentului 3. Ideea de bază este ceva de genul:

\

Cu alte cuvinte, \(\mathcal{A}\) duce la \(\mathcal{B}\), iar apoi \(\mathcal{B}\) duce la \(\mathcal{C}\). Cu modelele de mediere, presupunem o variabilă care intervine între covariabila normală \(\rightarrow\) și calea rezultatelor pe care am putea-o avea în cadrul regresiei standard, iar aceste modele ne permit să investigăm astfel de comportamente. În cazul de mai sus, variabila de intervenție, sau mediatorul, este \(\mathcal{B}\). Se întâmplă adesea să avem în continuare un efect direct al lui \(\mathcal{A}\) asupra lui \(\mathcal{C}\), dar, ca și în cazul modelului în general, acest lucru ar fi motivat teoretic.

Analiza de mediere este foarte populară în disciplinele din domeniul științelor sociale, deși nu se limitează în niciun caz la acestea, și se realizează de obicei sub masca modelării ecuațiilor structurale (SEM), care la rândul ei este o orientare specifică a modelelor grafice în general1. Modelul grafic al unui model de mediere ar putea arăta în felul următor.

În acest caz, a și b reflectă calea indirectă a efectului lui \(\mathrm{X}\) asupra rezultatului prin intermediul mediatorului, în timp ce c' este efectul direct al lui \(\mathrm{X}\ asupra rezultatului după ce calea indirectă a fost eliminată (c ar fi efectul înainte de a pune efectul indirect, iar c – c' este egal cu efectul indirect). Efectul total al \(\mathrm{X}\) este efectul indirect și direct combinat.

Ar trebui să notez câteva lucruri bazate pe ceea ce văd în consultanță în zeci de discipline. Pentru început, se pare că foarte puțini oameni care cred că au nevoie de un model de mediere chiar au nevoie. De exemplu, dacă nu vă puteți gândi la modelul dvs. în termeni temporali sau fizici, astfel încât \(\mathrm{X}\) să conducă în mod necesar la mediator, care apoi conduce în mod necesar la rezultat, probabil că nu aveți nevoie de un model de mediere. Dacă ați putea vedea săgețile mergând în ambele direcții, din nou, probabil că nu aveți nevoie de un astfel de model. De asemenea, dacă atunci când vă descrieți modelul, toată lumea crede că vorbiți despre o interacțiune (cunoscută și ca moderare), este posibil să nu aveți nevoie de acesta. Și, în cele din urmă, așa cum s-ar putea bănui, dacă nu există o corelație puternică între variabilele cheie (\(\mathrm{X}\)) și mediator (calea a) și dacă nu există o corelație puternică între mediator și rezultat(e) (calea b), probabil că nu aveți nevoie de acest lucru. Deși nimic nu vă va împiedica să faceți analiza de mediere, fără astfel de condiții prealabile, veți avea aproape sigur un model slab și probabil mai confuz decât ați fi avut altfel.

În concluzie, medierea funcționează cel mai bine atunci când există conexiuni cauzale puternic implicite între variabile. Chiar și atunci, un astfel de model ar trebui comparat cu un model mai simplu, fără mediere2. În orice caz, există câteva modalități foarte simple de a investiga astfel de modele în R, iar acesta este scopul aici, doar pentru a demonstra cum puteți începe.

Date

Pentru demonstrarea modelelor de mediere cu diferite pachete, vom folosi setul de date jobs care vine cu pachetul mediation. Iată descrierea.

Job Search Intervention Study (JOBS II). JOBS II este un experiment de teren randomizat care investighează eficacitatea unei intervenții de formare profesională asupra șomerilor. Programul este conceput nu numai pentru a crește rata de reangajare în rândul șomerilor, ci și pentru a îmbunătăți sănătatea mentală a persoanelor aflate în căutarea unui loc de muncă. În cadrul experimentului pe teren JOBS II, 1.801 șomeri au primit un chestionar de preselecție și au fost repartizați aleatoriu în grupuri de tratament și de control. Cei din grupul de tratament au participat la ateliere de pregătire profesională. În cadrul atelierelor, respondenții au învățat abilități de căutare a unui loc de muncă și strategii de adaptare pentru a face față eșecurilor în procesul de căutare a unui loc de muncă. Cei din grupul de control au primit o broșură care descria sfaturi de căutare a unui loc de muncă. În interviurile de urmărire, cele două variabile cheie de rezultat au fost o măsură continuă a simptomelor depresive pe baza listei de verificare a simptomelor Hopkins și o variabilă binară, reprezentând dacă respondentul a devenit angajat.

Iată o descriere a variabilelor din această demonstrație. Sunt disponibile și altele cu care ați putea dori să vă jucați.

  • econ_hard: Nivelul dificultăților economice înainte de tratament, cu valori de la 1 la 5.
  • sex: Variabilă indicatoare pentru sex. 1 = feminin
  • vârstă: Vârsta în ani.
  • educ: Factor cu cinci categorii pentru nivelul de educație.
  • job_seek: O scală continuă care măsoară nivelul de autoeficacitate în căutarea unui loc de muncă cu valori de la 1 la 5. Variabila mediatoare.
  • depress2: Măsură a simptomelor depresive post-tratament. Variabila de rezultat.
  • treat: Variabilă indicatoare pentru dacă participantul a fost selectat în mod aleatoriu pentru programul de formare JOBS II. 1 = atribuirea participării.
data(jobs, package = 'mediation')

Model

Date fiind aceste date, modelele pentru mediator și rezultat sunt următoarele:

\

Așa, ne așteptăm ca formarea competențelor profesionale să aibă un efect negativ asupra depresiei (adică o creștere a stării de bine), dar cel puțin o parte din aceasta s-ar datora unui efect pozitiv asupra căutării unui loc de muncă.

Ca un model grafic, l-am putea descrie succint după cum urmează.

Pachete

Vom examina următoarele pachete pentru a demonstra cum se poate efectua o analiză de mediere în R:

  • mediation
  • lavaan
  • psych
  • brms

În timp ce acestea vor fi în centrul atenției, voi nota și alte câteva alternative, inclusiv Python și Stata.

mediation

Vom începe cu pachetul de mediere, deoarece, practic, nu necesită mai multe abilități de programare pentru a conduce decât cele pe care cineva le posedă deja din rularea modelelor de regresie standard în R. Pachetul oferă efectul mediu de mediere cauzală, definit după cum urmează din fișierul de ajutor și din articolele lui Imai3:

Efectul mediu de mediere cauzală (ACME) reprezintă diferența așteptată în rezultatul potențial atunci când mediatorul a luat valoarea care s-ar realiza în condiția de tratament față de condiția de control, în timp ce statutul de tratament în sine este menținut constant.

Rețineți cum această definiție se concentrează pe valorile așteptate sau prezise condiționate de valoarea tratamentului. Această noțiune de contrafactualitate, sau cum ar arăta observația în condiții opuse, are o istorie îndelungată în modelare în acest punct. Gândiți-vă în felul următor: dacă cineva face parte din grupul de tratament, ar avea o anumită valoare pentru mediator și, având în vedere acest lucru, ar avea apoi o anumită valoare așteptată pentru rezultat. Cu toate acestea, am putea să presupunem că aceeași observație se află și în grupul de control și să evaluăm efectul asupra rezultatului prin intermediul mediatorului la fel de bine. Putem evalua rezultatele potențiale, menținând constant tratamentul. Gândindu-ne la modificările rezultatelor având în vedere valoarea mediatorului nu face nicio presupunere cu privire la tipul de model. Acesta este modul în care pachetul de mediere este capabil să încorporeze diferite modele pentru mediator vs. rezultat. De exemplu, mediatorul ar putea fi binar, necesitând un model de regresie logistică, în timp ce modelul rezultatului ar putea fi un model de supraviețuire.

În exemplul nostru, vom rămâne la modele liniare standard (normale). Rețineți, de asemenea, că, deși tratamentul nostru este o variabilă binară, acest lucru se generalizează la cazul continuu, în care luăm în considerare rezultatul unei mișcări de o unitate pe „tratament”. Pentru ca pachetul de mediere să funcționeze, pur și simplu rulăm modelele noastre respective pentru mediator și rezultat, apoi folosim funcția mediate pentru a obține rezultatul final.

library(mediation)model_mediator <- lm(job_seek ~ treat + econ_hard + sex + age, data = jobs)model_outcome <- lm(depress2 ~ treat + econ_hard + sex + age + job_seek, data = jobs)# Estimation via quasi-Bayesian approximation?mediatemediation_result <- mediate( model_mediator, model_outcome, sims = 500, treat = "treat", mediator = "job_seek")detach(package:mediation)detach(package:MASS)
summary(mediation_result)plot(mediation_result)
Estimate 95% CI Lower 95% CI Upper p-.value
ACME -0.016 -0,038 0,009 0,220
ADE -0.045 -0,127 0,047 0,292
Efect total -0.061 -0,149 0,027 0,188
Prop. mediată 0.226 -3,222 1,596 0,344

Rezultatele de mai sus demonstrează că ACME nu se deosebește statistic de zero, sau de lipsa medierii. Efectul direct mediu este negativ, dar, de asemenea, nu este notabil din punct de vedere statistic și nici efectul total (efect indirect + efect direct). De asemenea, este furnizat soi disant „proporția mediată”, care este raportul dintre efectul indirect și cel total. Cu toate acestea, aceasta nu este o proporție și poate fi chiar negativă, deci este în mare parte un număr lipsit de sens.

Pros

  • Modeluri și sintaxă standard R
  • Multe tipuri de modele atât pentru mediator, cât și pentru rezultat
  • Furnizează mai multe rezultate simultan
  • Documentație bună și articole asociate sunt disponibile gratuit
  • Poate face mediere „moderată”

Limitări

    .

  • Utilizarea MASS4
  • Modeluri simple cu efecte aleatorii
  • Funcționalitatea poate fi limitată cu unele complexități ale modelului
  • Nu există capacități de variabile latente

lavaan

În cazul specific în care ambele modele de mediere și de rezultat sunt modele liniare standard cu o distribuție normală pentru variabila țintă, efectul indirect este echivalent cu produsul dintre traseele a și b din diagrama anterioară. Efectul direct este calea c'. O comparație a efectului direct de sine stătător, pe care l-am putea numi c, față de acest efect direct estimat în modelul de mediere c', este de așa natură încât c - c' = a*b. Ceea ce a fost menționat mai devreme ar putea fi acum mai clar, dacă fie a, fie b sunt aproape zero, atunci efectul indirect nu poate fi decât aproape zero, astfel încât este prudent să se investigheze în prealabil astfel de relații.

Această abordare a produsului de căi (sau a diferenței de coeficienți) este cea pe care o vom adopta cu pachetul lavaan și, de fapt, la momentul scrierii acestui articol, aceasta este singura noastră modalitate de a proceda. lavaan este orientat în mod special către modelarea ecuațiilor structurale, cum ar fi analiza factorială, modelele de creștere și modelele de mediere, cum este cel pe care îl realizăm aici, și este foarte recomandat pentru astfel de modele. Deși este limitat la cazul modelului liniar standard pentru a evalua medierea, este singurul dintre instrumentele noastre care poate încorpora cu ușurință variabile latente5. De exemplu, am putea avea rezultatul nostru privind depresia ca o variabilă latentă care stă la baza itemilor chestionarului individual. În plus, am putea încorpora, de asemenea, mai mulți mediatori și mai multe rezultate.

Pentru a păstra lucrurile așa cum am discutat, voi eticheta căile a, b și c' în lavaan în funcție de modul în care au fost reprezentate anterior. În rest, lavaan este foarte ușor de utilizat și, în cazul variabilelor observate, utilizează notația standard a formulelor R pentru modele. Dincolo de aceasta, definim efectele de interes pe care dorim să le calculăm cu operatorul :=. Specificăm modelul în întregime ca un simplu șir de caractere, apoi folosim funcția sem pentru a face analiza. 

library(lavaan)sem_model = ' job_seek ~ a*treat + econ_hard + sex + age depress2 ~ c*treat + econ_hard + sex + age + b*job_seek # direct effect direct := c # indirect effect indirect := a*b # total effect total := c + (a*b)'model_sem = sem(sem_model, data=jobs, se='boot', bootstrap=500)summary(model_sem, rsq=T) # compare with ACME in mediation
lavaan 0.6-6 ended normally after 25 iterations Estimator ML Optimization method NLMINB Number of free parameters 11 Number of observations 899 Model Test User Model: Test statistic 0.000 Degrees of freedom 0Parameter Estimates: Standard errors Bootstrap Number of requested bootstrap draws 500 Number of successful bootstrap draws 500Regressions: Estimate Std.Err z-value P(>|z|) job_seek ~ treat (a) 0.066 0.049 1.332 0.183 econ_hard 0.053 0.024 2.242 0.025 sex -0.008 0.047 -0.163 0.870 age 0.005 0.002 1.934 0.053 depress2 ~ treat (c) -0.040 0.044 -0.905 0.365 econ_hard 0.149 0.022 6.908 0.000 sex 0.107 0.038 2.831 0.005 age 0.001 0.002 0.332 0.740 job_seek (b) -0.240 0.030 -8.079 0.000Variances: Estimate Std.Err z-value P(>|z|) .job_seek 0.524 0.030 17.610 0.000 .depress2 0.373 0.022 17.178 0.000R-Square: Estimate job_seek 0.011 depress2 0.120Defined Parameters: Estimate Std.Err z-value P(>|z|) direct -0.040 0.045 -0.904 0.366 indirect -0.016 0.012 -1.324 0.185 total -0.056 0.046 -1.224 0.221

Vezi aceeași ieșire anterioară și putem compara parametrul nostru indirect cu ACME pe care îl aveam înainte, efectul direct este comparat cu ADE, iar total se compară cu efectul total anterior. Valorile sunt în esență aceleași.

Rețineți, de asemenea, că ieșirea arată valoarea \(R^2\) pentru ambele modele. În cazul job_seek, putem vedea că motivul pentru care nu găsim prea multe în ceea ce privește medierea este acela că, pentru început, covariatele implicate nu explică nicio variație a mediatorului. O investigație preliminară ne-ar fi scutit de probleme în acest caz.

Pros

  • Puteți gestiona mai mulți mediatori
  • Puteți gestiona mai multe „tratamente”
  • Puteți gestiona mai multe rezultate
  • Puteți utiliza variabile latente
  • Un oarecare sprijin pe mai multe niveluri
  • Puteți face mediere moderată și mediată moderare (deși nu pentru variabile latente)

Limitări

  • Este necesară o codificare suplimentară pentru a estima efectul indirect
  • Efecte aleatorii unice
  • În timp ce modelele ar putea încorpora variabile binare sau ordinale pentru mediator/rezultate, nu există o modalitate directă de a calcula efectul indirect în maniera pachetului de mediere în aceste contexte.

piecewiseSEM

Pachetul piecewiseSEM funcționează foarte asemănător cu pachetul de mediere. Ceea ce este interesant în raport cu pachetul de mediere este că piecewiseSEM poate gestiona tipuri suplimentare de modele, precum și să ofere ieșiri suplimentare (de exemplu, rezultate standardizate), opțiuni suplimentare (de exemplu, multigrup, reziduuri corelate) și vizualizarea modelului.

library(piecewiseSEM)model_mediator <- lm(job_seek ~ treat + econ_hard + sex + age, data = jobs)model_outcome <- lm(depress2 ~ treat + econ_hard + sex + age + job_seek, data = jobs)mediation_result <- psem(model_mediator, model_outcome, data = jobs)summary(mediation_result)
Structural Equation Model of mediation_result Call: job_seek ~ treat + econ_hard + sex + age depress2 ~ treat + econ_hard + sex + age + job_seek AIC BIC 26.000 88.417---Tests of directed separation: No independence claims present. Tests of directed separation not possible.Global goodness-of-fit: Fisher's C = 0 with P-value = 1 and on 0 degrees of freedom---Coefficients: Response Predictor Estimate Std.Error DF Crit.Value P.Value Std.Estimate job_seek treat 0.0656 0.0515 894 1.2748 0.2027 0.0425 job_seek econ_hard 0.0532 0.0246 894 2.1612 0.0309 0.0720 * job_seek sex -0.0076 0.0487 894 -0.1567 0.8755 -0.0052 job_seek age 0.0046 0.0023 894 1.9779 0.0482 0.0658 * depress2 treat -0.0403 0.0435 893 -0.9255 0.3550 -0.0291 depress2 econ_hard 0.1485 0.0208 893 7.1323 0.0000 0.2248 *** depress2 sex 0.1068 0.0411 893 2.5957 0.0096 0.0818 ** depress2 age 0.0006 0.0020 893 0.3306 0.7410 0.0104 depress2 job_seek -0.2400 0.0282 893 -8.4960 0.0000 -0.2682 *** Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05---Individual R-squared: Response method R.squared job_seek none 0.01 depress2 none 0.12

Ne putem folosi de capacitățile sale de reprezentare grafică pentru a crea o vizualizare rapidă a modelului.

plot(mediation_result)

Din păcate, în prezent nu există o modalitate automată de a calcula efectele indirecte, astfel încât ar trebui să se bootstrap rezultatele manual.

Pros

  • Modeluri și sintaxă standard R
  • Multe tipuri de modele atât pentru mediator cât și pentru rezultat
  • Câteva rezultate de tip SEM (de ex. potrivire, coeficienți standardizați, AIC)
  • Traiect rapid al rezultatelor
  • Poate gestiona mai mulți mediatori, „tratamente”, și rezultate

Limitări

  • Nu calculează automat efectele indirecte
  • Nu dispune de capacități de variabile latente

psihologie

Pachetul psych profită de faptul că, în cazul modelului liniar standard, se pot obține rezultatele prin intermediul modelelor de regresie adecvate bazate doar pe matricile de covarianță. Este foarte asemănător cu lavaan, deși folosește o abordare prin metoda celor mai mici pătrate ordinare, spre deosebire de cea de maximă verosimilitate. Ceea ce este frumos aici este o sintaxă care vă permite să vă concentrați doar asupra efectului de interes sau să includeți totul, ceea ce este frumos dacă ați fi interesat și de efectele indirecte pentru dificultățile economice, vârstă și sex.

Pentru acest demo vom folosi versiunea curățată folosind -, în loc de +, pentru efectele non-tratament. Acest lucru înseamnă doar că acestea sunt incluse cu modelele, dar rezultatele nu sunt afișate în ceea ce le privește. Mediatorul este identificat cu (). Un alt bonus este un grafic rapid al rezultatelor, care arată diferența dintre efectele directe neajustate și ajustate, precum și intervalul bootstrap corespunzător. 

library(psych)mediation_psych = mediate( depress2 ~ treat + (job_seek) - econ_hard - sex - age, data = jobs, n.iter = 500)

mediation_psych
Mediation/Moderation Analysis Call: mediate(y = depress2 ~ treat + (job_seek) - econ_hard - sex - age, data = jobs, n.iter = 500)The DV (Y) was depress2* . The IV (X) was treat* . The mediating variable(s) = job_seek* . Variable(s) partialled out were econ_hard sex ageTotal effect(c) of treat* on depress2* = -0.06 S.E. = 0.05 t = -1.24 df= 895 with p = 0.21Direct effect (c') of treat* on depress2* removing job_seek* = -0.04 S.E. = 0.15 t = 14.91 df= 893 with p = 4.6e-45Indirect effect (ab) of treat* on depress2* through job_seek* = -0.02 Mean bootstrapped indirect effect = -0.02 with standard error = 0.01 Lower CI = -0.04 Upper CI = 0.01R = 1.07 R2 = 1.15 F = -3510.4 on 2 and 893 DF p-value: 1 To see the longer output, specify short = FALSE in the print statement or ask for the summary
summary(mediation_psych)
Call: mediate(y = depress2 ~ treat + (job_seek) - econ_hard - sex - age, data = jobs, n.iter = 500)Direct effect estimates (traditional regression) (c') depress2* se t df ProbIntercept 2.21 0.15 14.91 893 4.60e-45treat -0.04 0.04 -0.93 893 3.55e-01job_seek -0.24 0.03 -8.50 893 8.14e-17R = 1.07 R2 = 1.15 F = -3510.4 on 2 and 893 DF p-value: 1 Total effect estimates (c) depress2* se t df Probtreat -0.06 0.05 -1.24 895 0.215 'a' effect estimates job_seek se t df ProbIntercept 3.67 0.13 29.33 894 5.65e-133treat 0.07 0.05 1.27 894 2.03e-01 'b' effect estimates depress2* se t df Probjob_seek -0.24 0.03 -8.5 894 7.83e-17 'ab' effect estimates (through mediators) depress2* boot sd lower uppertreat -0.02 -0.02 0.01 -0.04 0.01

Aceleași rezultate, un ambalaj diferit, dar probabil cea mai ușoară cale de până acum, deoarece a necesitat doar un apel de funcție. Pachetul psihologic gestionează, de asemenea, mediatori și rezultate multiple, ca un bonus.

Pros

  • Cea mai simplă sintaxă, practic un model pe o singură linie
  • Practicare rapidă a rezultatelor
  • Poate gestiona mai mulți mediatori, „tratamente”, și rezultate
  • Poate face mediere „moderată”

Limitări

  • Limitat la modelul liniar standard (lm)
  • Utilizarea MASS

brms

Pentru următoarea demonstrație, ajungem la ceea ce consider a fi cel mai puternic pachet, brms. Numele vine de la Bayesian Regression Modeling with Stan, iar Stan este un puternic limbaj de programare probabilistică pentru analiza bayesiană. Nu voi intra în detalii despre analiza bayesiană, dar nu ezitați să consultați documentul meu care o face.

În general, facem ca și până acum, specificând modelul mediatorului și modelul rezultatului. brms nu face nimic special pentru analiza de mediere, dar funcția sa de ipoteză ne poate permite să testăm abordarea produsului traseelor. Mai mult decât atât, pachetul sjstats va furniza în esență rezultatele în același mod în care pachetul de mediere o face pentru noi și, de altfel, pachetul de mediere este practic o încercare de a găsi o soluție bayesiană folosind oricum metode frequentiste. Dacă am avea distribuții diferite pentru rezultat și mediator, ne-ar fi relativ ușor să obținem aceste valori medii de predicție și diferențele lor, deoarece abordările bayesiene se gândesc întotdeauna la distribuțiile predictive posterioare. În orice caz, iată codul.

library(brms)model_mediator <- bf(job_seek ~ treat + econ_hard + sex + age)model_outcome <- bf(depress2 ~ treat + job_seek + econ_hard + sex + age)med_result = brm( model_mediator + model_outcome + set_rescor(FALSE), data = jobs)save(med_result, file = 'data/mediation_brms.RData')
load('data/mediation_brms.RData')summary(med_result)
 Family: MV(gaussian, gaussian) Links: mu = identity; sigma = identity mu = identity; sigma = identity Formula: job_seek ~ treat + econ_hard + sex + age depress2 ~ treat + job_seek + econ_hard + sex + age Data: jobs (Number of observations: 899) Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1; total post-warmup samples = 4000Population-Level Effects: Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESSjobseek_Intercept 3.67 0.12 3.43 3.91 1.00 6699 3749depress2_Intercept 2.21 0.15 1.92 2.50 1.00 6174 3091jobseek_treat 0.07 0.05 -0.03 0.17 1.00 6322 2709jobseek_econ_hard 0.05 0.02 0.00 0.10 1.00 6266 2656jobseek_sex -0.01 0.05 -0.10 0.09 1.00 5741 2655jobseek_age 0.00 0.00 0.00 0.01 1.00 6539 2846depress2_treat -0.04 0.04 -0.12 0.04 1.00 5458 3102depress2_job_seek -0.24 0.03 -0.30 -0.18 1.00 5950 2938depress2_econ_hard 0.15 0.02 0.11 0.19 1.00 7543 3102depress2_sex 0.11 0.04 0.03 0.19 1.00 5599 2699depress2_age 0.00 0.00 -0.00 0.00 1.00 4555 2887Family Specific Parameters: Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESSsigma_jobseek 0.73 0.02 0.69 0.76 1.00 6639 3276sigma_depress2 0.61 0.01 0.59 0.64 1.00 6145 2987Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESSand Tail_ESS are effective sample size measures, and Rhat is the potentialscale reduction factor on split chains (at convergence, Rhat = 1).
# using brms we can calculate the indirect effect as follows# hypothesis(med_result, 'jobseek_treat*depress2_job_seek = 0')# sjstats provides similar printing as the mediation package# print(sjstats::mediation(med_result), digits=4)sjstats::mediation(med_result) %>% kable_df()
efect valoare hdi.low hdi.high
direct -0.039 -0,112 0,031
indirect -0,015 -0,036 0,005
mediator -0.240 -0,286 -0,193
total -0,055 -0,133 0,017
proporția mediată 0.277 -0,813 1,366

În ieșire, tot ceea ce are jobseek_* este un rezultat pentru modelul mediator, în timp ce depress2_* este pentru rezultat. Avem aceeași poveste veche în acest punct, dar cu abordarea bayesiană avem mai multe lucruri amuzante la care să ne uităm. De exemplu, putem vedea că, de fapt, nu capturăm bine asimetria rezultatului depresiei. Valorile noastre previzionate în comparație cu cele observate nu prea se potrivesc. Suntem un pic mai bine pentru mediator, dar poate că suntem încă un pic cam mari cu unele dintre predicțiile noastre bazate pe model.

pp_check(med_result, resp = 'depress2') + ggtitle('Depression Outcome')
pp_check(med_result, resp = 'depress2') + ggtitle('Depression Outcome')

pp_check(med_result, resp = 'jobseek') + ggtitle('Mediator')

Pros

  • Sintaxă directă
  • Extrem de puternică… Modelele sunt în mare parte limitate de imaginația fiecăruia
  • În principiu, face ceea ce pachetul de mediere aproximează
  • Toate avantajele inferenței bayesiene: diagnostice, verificări predictive posterioare, compararea modelelor, etc.

Limitări

  • Mai ușor de estimat
  • Calculele „de mână” necesare pentru a depăși modelul liniar standard, dar aceasta este deja o abordare comună din perspectiva bayesiană
  • Este necesar un anumit confort cu abordarea bayesiană

Mai multă complexitate

Câteva dintre pachetele menționate pot gestiona modele mai complexe sau pot oferi abordări suplimentare pentru a investiga efectele indirecte.

Interacțiuni

Câteva modele implică interacțiuni fie pentru modelul de mediere, fie pentru rezultat și, din nefericire, acest lucru este adesea denumit moderare mediată sau mediere moderată. Personal, nu văd avantajul de a da denumiri ambigue la ceea ce, altfel, ar putea fi un concept simplu (chiar dacă încă nu este un model atât de simplu), dar această navă a navigat de mult timp. Nu am de gând să intru în detalii, dar ideea este că s-ar putea să aveți un termen de interacțiune undeva în model, iar interacțiunea ar putea implica variabila de tratament, mediatorul sau ambele.

Suficient să spunem că, din moment ce folosim instrumente standard de modelare precum lm și extensii ale acestuia, încorporarea interacțiunilor este trivială pentru toate pachetele de mai sus, dar abordarea de tip produs de trasee nu este valabilă (a*b != c').

Modeluri liniare generalizate

În unele cazuri mediatorul sau rezultatul nostru poate fi binar, numărat, sau ceva în care presupunerea unei distribuții normale ar putea să nu fie cea mai bună idee. Sau am putea dori să investigăm relații neliniare între tratament/mediator/rezultat. Sau am putea avea date care au observații corelate, cum ar fi măsurători repetate sau altele similare. Pachetul de mediere se mândrește în special cu acest lucru, dar brms poate face tot ceea ce poate face și chiar mai mult, deși s-ar putea să trebuiască să depuneți ceva mai multă muncă pentru a calcula efectiv rezultatul. lavaan poate de fapt să facă un set limitat de modele pentru variabile binare și ordinale, dar obținerea estimării indirecte adecvate ar necesita o abordare manuală foarte plictisitoare.

Date lipsă

De multe ori, atunci când avem de-a face cu astfel de date, în special în științele sociale, datele lipsesc adesea pentru oricare dintre covariate. Uneori putem renunța la acestea dacă nu sunt prea multe, dar în alte cazuri vom dori să facem ceva în acest sens. Pachetele lavaan, psych și brms oferă una sau mai multe modalități de a face față situației (de exemplu, imputare multiplă).

Alternative

Am descris modelele ca rețele de noduri, cu arce/izvoare/căi care le conectează. Discuția noastră se învârte în jurul a ceea ce se numește grafuri aciclice dirijate (DAG) în care săgețile pot merge doar într-o singură direcție, fără bucle de feedback. Rezultatul oricărei variabile de rezultat este o funcție a săgeților care o precedă și este condiționat independent de celelalte. Unele modele teoretice pot relaxa acest lucru, iar altele pot să nu aibă deloc săgeți, adică sunt nedirecționate, astfel încât suntem interesați doar de conexiuni (de exemplu, în cazul unor rețele sociale).

bnlearn

Pachetul bnlearn permite investigarea grafurilor direcționate, parțial direcționate și nedirecționate. În ceea ce privește DAG-urile, îl putem folosi pentru a duplica, în esență, modelele de mediere pe care le-am discutat. Totuși, partea bună este că acest pachet va testa în mod eficient căile de includere în loc să le presupună, dar putem impune în continuare constrângeri teoretice, după cum este necesar. Nu numai că putem căuta apoi căile de interes într-un mod principial cu ajutorul rețelelor bayesiene și al teoriei grafurilor cauzale a lui Pearl ca bază, dar vom avea, de asemenea, instrumente pentru a evita în continuare supraadaptarea prin intermediul validării încrucișate.

Pentru modelul inițial, ne vom asigura că există căi între tratament – mediator, tratament – rezultat și mediator – rezultat (lista albă). Vom refuza căile fără sens, cum ar fi faptul de a avea săgeți către tratament (care a fost alocat aleatoriu), sexul, dificultățile economice și vârsta (lista neagră). În rest, vom vedea ce sugerează datele. 

whitelist = data.frame( from = c('treat', 'treat', 'job_seek'), to = c('job_seek', 'depress2', 'depress2'))blacklist = expand.grid( from = colnames(mediation_result$model.y$model), to = c('treat', 'sex', 'age', 'econ_hard'))# For simpler output we'll use treatment and sex as numeric (explained later)library(dplyr)jobs_trim = jobs %>% select(depress2, treat, econ_hard, sex, age, job_seek) %>% mutate( treat = as.numeric(jobs$treat), sex = as.numeric(jobs$sex) )# extract path coefficients if desired# parameters = bn.fit(model, jobs_trim)# parameters$job_seek# parameters$econ_hard# parameters$depress2
library(bnlearn)model = gs(jobs_trim, whitelist = whitelist, blacklist = blacklist)plot(model)

Vezi în grafic că lucrurile s-au schimbat puțin. De exemplu, vârsta se leagă acum doar de autoeficacitatea în căutarea unui loc de muncă, iar sexul are efect doar asupra depresiei.

Dacă restrângem căile să fie doar ceea ce sunt în exemplele noastre anterioare, vom obține aceleași rezultate.

library(bnlearn)whitelist = data.frame( from = c('treat', 'age', 'sex', 'econ_hard', 'treat', 'job_seek', 'age', 'sex', 'econ_hard'), to = c('job_seek', 'job_seek','job_seek','job_seek', 'depress2', 'depress2', 'depress2', 'depress2', 'depress2'))blacklist = expand.grid( from = colnames(mediation_result$model.y$model), to = c('treat', 'sex', 'age', 'econ_hard')) model = gs(jobs_trim, whitelist = whitelist, blacklist = blacklist)plot(model)

parameters = bn.fit(model, jobs_trim)parameters$depress2$coefficients
 (Intercept) treat econ_hard sex age job_seek 2.2076414333 -0.0402647000 0.1485433818 0.1068048699 0.0006488642 -0.2399549527 
parameters$job_seek$coefficients
 (Intercept) treat econ_hard sex age 3.670584908 0.065615003 0.053162413 -0.007637336 0.004586492

Principalul lucru care trebuie remarcat este că parametrii estimați sunt egali cu ceea ce am obținut cu pachetele anterioare. În esență, este echivalent cu utilizarea lavaan cu estimatorul implicit de maximă verosimilitate.

Dacă folosim tratamentul și sexul ca factori, bnlearn va produce modele condiționate care sunt diferite în funcție de valoarea factorului luată. Cu alte cuvinte, se va avea un model separat pentru când treatment == 'treatment' și unul pentru când treatment == control. În cazul nostru, acest lucru ar fi identic cu a permite ca totul să interacționeze cu tratamentul, de exemplu lm( job_seek ~ treat * (econ_hard + sex + age)), și la fel pentru modelul de depresie. Acest lucru s-ar extinde la potențial orice variabilă binară (de exemplu, inclusiv sexul). Dacă mediatorul este o variabilă binară, acest lucru este probabil ceea ce am dori să facem.

Python

Directorul CCSAR Kerby Shedden a ținut un atelier Python privind modelele de mediere, așa că arăt implementarea statsmodels aici. Aceasta urmează abordarea Imai și astfel poate fi văzută ca versiunea Python a pachetului de mediere. Rezultatul este în esență același cu cel pe care l-ați avea folosind tratamentul ca variabilă factorială, unde obțineți rezultate separate pentru fiecare categorie de tratament. Acest lucru nu este necesar pentru demonstrația noastră, așa că puteți compara pur și simplu rezultatele „medii” cu rezultatele anterioare ale pachetului de mediere. 

import statsmodels.api as smfrom statsmodels.stats.mediation import Mediationimport numpy as npimport pandas as pdoutcome_model = sm.OLS.from_formula("depress2 ~ treat + econ_hard + sex + age + job_seek", data = jobs)mediator_model = sm.OLS.from_formula("job_seek ~ treat + econ_hard + sex + age", data = jobs)med = Mediation(outcome_model, mediator_model, "treat", "job_seek")med_result = med.fit(n_rep = 500)print(np.round(med_result.summary(), decimals = 3))
 Estimate Lower CI bound Upper CI bound P-valueACME (control) -0.016 -0.048 0.014 0.332ACME (treated) -0.016 -0.048 0.014 0.332ADE (control) -0.043 -0.130 0.044 0.308ADE (treated) -0.043 -0.130 0.044 0.308Total effect -0.059 -0.144 0.029 0.208Prop. mediated (control) 0.241 -1.710 2.254 0.364Prop. mediated (treated) 0.241 -1.710 2.254 0.364ACME (average) -0.016 -0.048 0.014 0.332ADE (average) -0.043 -0.130 0.044 0.308Prop. mediated (average) 0.241 -1.710 2.254 0.364

Stata

În cele din urmă, furnizez o opțiune în Stata, folosind comanda sem a acestuia. Stata facilitează obținerea efectelor indirecte în acest exemplu, dar face acest lucru pentru fiecare covariantă, astfel încât rezultatul este cel puțin un pic verbos6. Pentru cei care lucrează cu Stata, nu au nevoie de un pachet SEM separat pentru a obține acest tip de rezultate.

use "data\jobs.dta"sem (job_seek <- treat econ_hard sex age) (depress2 <- treat econ_hard sex age job_seek), cformat(%9.3f) pformat(%5.2f)estat teffects, compact cformat(%9.3f) pformat(%5.2f)

Rezumat

Modelurile cu efecte indirecte necesită o analiză teoretică atentă pentru a fi utilizate pentru analiza datelor. Cu toate acestea, dacă modelul este adecvat pentru situația datelor dumneavoastră, este destul de ușor să obțineți rezultate dintr-o varietate de pachete din R. Mai mult, nu este necesar să se utilizeze un pachet de modelare a ecuațiilor structurale pentru a efectua o analiză cu efecte indirecte și, de fapt, se poate ajunge departe folosind sintaxa standard R. Pentru variabilele strict observate, adică fără variabile latente, nu este necesar sau chiar recomandat niciun instrument SEM.

\

Compararea pachetelor rezumată

Următorul tabel poate ajuta pe cineva să decidă ce pachet să folosească pentru nevoile sale, având în vedere considerațiile sale teoretice.

.

.

.

.

.

.

mediation lavaan piecewiseSEM psych brms
Automat – – -*
Tratamente multiple☺ – – – – – – – –
Multiplu Mediatori – – – – – – –
Rezultate multiple – – – – –
Beyond SLM† – – – – – –
Efecte aleatorii – – – – – – –
Valori lipsă – -* – –
Variabile latente – -*
* aproximativ, cu unele avertismente
☺ Poate necesita reluarea unor aspecte ale modelului
† Model liniar standard, așa cum a fost estimat de lm

  1. Am un document mult mai detaliat despre SEM, inclusiv despre analiza de mediere.︎

  2. Dintr-un motiv oarecare nu prea vezi acest lucru în practică și ne întrebăm ce s-a făcut pentru ca datele să se preteze la un astfel de model dacă acesta nu era justificat.︎

  3. Imai își pune la dispoziție articolele pe site-ul său.︎

  4. MASS a fost înlocuit de alții de mai bine de un deceniu în acest moment, iar în cea mai mare parte tinde doar să vă strice tidyverse și alte pachete atunci când este încărcat. Este un pachet bun (și a fost grozav pe vremuri), dar dacă vreți să îl folosiți într-un pachet, ar fi bine să nu îl încărcați (sau alte pachete) în mediu doar pentru a folosi o funcție sau două. De cele mai multe ori îl văd folosit doar pentru mvrnorm (distribuție normală multivariată) și glm.nb, dar există alte pachete cu această funcționalitate care ar oferi beneficii suplimentare, și nu funcții dplyr mascate, care sunt printre cele mai utilizate în comunitatea R.︎

  5. brms lucrează la asta.︎

  6. Opțiunile din cod sunt acolo pentru a suprima/minimiza ceea ce poate fi.︎

Share:

Lasă un răspuns

Adresa ta de email nu va fi publicată.