Gaat opmerking in Bain tutorial ook voor ons op??
Ik hoop het niet...

opmerking in Bain tutorial:
Two-way or higher order ANOVA's can only be handled by recoding all factors into one factor. If, or example, there is a factor sex with levels man and woman, and a factor age with levels young and old, these have to be recoded in a new factor sexage with levels manyoung, manold, womanyoung, womanold

Ik wil bijv, als lik = max lik is, iets aangeven.

Ik heb dit als issue, zodat ik kan zien dat je dit gedaan hebt.
Dan weet ik wanneer ik met de tutorials bezig kan (omdat ik dan die van Herbert kan nabootsen).

De gorica package maakt in de comp() functie grotendeels gebruik van de complement code uit restriktor zonder hier een acknowledgement te plaatsen. Het is wel netjes om dit aan het pakket toe te voegen, bijvoorbeeld als comment. Daarnaast, staat in de code nog verwijzingen naar restriktor:

stop("Restriktor ERROR: no complement exists for the unrestricted hypothesis.")
en
stop("Restriktor ERROR: you might have found a bug, please contact me!")

Zou je dit willen aanpassen.

Uitbreidingen die ik graag nog zie:

• Include goric values (voor nu iig wanneer Multivariaat normal linear models – tzt wellicht voor meer type modellen, bijv ook SEM).

4. Kan jij dat doen obv de goric R package?

In at least one of the testthat examples (test-sem_example.R), the following lavaan syntax is used:

model <- '
Cry = ~ y1 + y2 + y3 + y4
Fld = ~ y2 + y3 + y5 + y6 + y7 + y8
Cry ~ edc + age
Fld ~ edc + age
'

Note that in the first two syntax lines, there is space beteen the '=' and '~' symbols. However, lavaan's new parser will no longer allow for this. Is it possible to remove the space in this (and perhaps other) example(s)? Otherwise, it will break the testthat tests and this will prevent the new version of lavaan (0.6-17) to get through the CRAN checks.

De gorica functie geeft misfit en complexity weer. Op zich is het idee prima, maar ik zou toch graag het volgende zien: fit = loglik en complexity = penalty.
Heel evt kan er in een note ergens dat gorica = -2 * fit + 2 * complexity of gorica = -2 * loglik + 2 * penalty.

• Straks contingency table code (levert Yasin aan) – dus code specifiek voor bepaalde type modellen (contingency tables en evt multinomial regression).
  Btw, then link functions (e.g., from prob’s to odds ratios), but can be easily handed in text.

Also incorporate the functionality of evaluating a hypothesis against its complement (and if more than one substantive hypothesis: also the complement of the set).

Ik kon zo snel niet zien dat er een beschrijving staat bij de gorica mbt het lavaan object. Zo weet de gebruiker ook niet dat er een 'standardize=TRUE' optie/command is...
Zie evt code hieronder.

Load the gorica and lavaan libraries.

if (!require("gorica")) install.packages("gorica") # install this package first (once)
if (!require("lavaan")) install.packages("lavaan") # install this package first (once)
library(gorica)
library(lavaan)

Specify the latent regression model

model2 <- '
A =~ Ab + Al + Af + An + Ar + Ac
B =~ Bb + Bl + Bf + Bn + Br + Bc
A ~ B + age + peabody
'

Fit the latent regression model using the lavaan sem function

fit2 <- sem(model2, data = sesamesim, std.lv = TRUE)

Formulate hypotheses

hypotheses2 <- "
AB > Apeabody = Aage = 0;
AB > Apeabody > Aage = 0;
AB > Apeabody > A~age > 0
"

Call gorica

Note: because more than 1 hypothesis, we cannot use: comparison = "complement"

set.seed(100)
out2_c <- gorica(fit2, hypotheses2, standardize=TRUE)
out2_c

Op dit moment is er alleen code om de resultaten van 1 hypothese vs zijn complement te bepalen en niet voor een set van hypotheses. Toen ik het deed voor meer dan een hypothese, dan geeft het wel output. Op zich geeft het de output met de unconstrained, alleen is er geen melding dat dit gebeurd en daarnaast staat er in de output ook dat het de complement is...
Kan je een foutmelding inbouwen. Als je daarnaast ook de resultaten wilt geven met de unconstrained dan moet de hypothese daar aangepast worden (niet Hc en ook niet de omschrijving - zie hieronder).

Hieronder vind je code en output.

Load the gorica and lavaan libraries.

if (!require("gorica")) install.packages("gorica") # install this package first (once)
if (!require("lavaan")) install.packages("lavaan") # install this package first (once)
library(gorica)
library(lavaan)

Specify the latent regression model

model2 <- '
A =~ Ab + Al + Af + An + Ar + Ac
B =~ Bb + Bl + Bf + Bn + Br + Bc
A ~ B + age + peabody
'

Fit the latent regression model using the lavaan sem function

fit2 <- sem(model2, data = sesamesim, std.lv = TRUE)

Formulate hypotheses

hypotheses2 <- "
AB > Apeabody = Aage = 0;
AB > Apeabody > Aage = 0;
AB > Apeabody > A~age > 0
"

Call gorica

Note: because more than 1 hypothesis, we cannot use: comparison = "complement"

set.seed(100)
out2_c <- gorica(fit2, hypotheses2, comparison = "complement", standardize=TRUE)
out2_c

OUTPUT:

Informative hypothesis test for an object of class lavaan:

loglik penalty gorica gorica_weights
H1 6.836 0.504 -12.663 0.370
H2 6.836 0.691 -12.290 0.307
H3 6.836 0.823 -12.026 0.269
Hc 6.894 2.496 -8.798 0.054

Hypotheses:
H1: AB>Apeabody=Aage=0
H2: AB>Apeabody>Aage=0
H3: AB>Apeabody>A~age>0
Hc: Complement of the hypothesis

How to handle this? lm includes data, which can be standardized and the model re-run. Other model objects do not include data.

Rebecca: Include ‘c’ as option, then range restristrictions (like 0 < mu1 < 5) possible (but dependent on c). Kan ik tzt meer uitleg over geven. Code hiervoor ligt er nog niet, dit is voor in de verdere toekomst 😊.

• You also make this available for Mplus input (via Mplus automation) and make sure the estimates have names assigned etc.
  Could you also then standardize the data and re-run the model?
• You also make this available for Lavaan objects (like you will do for Bain).

6. De hypotheses kunnen nu in termen van namen worden gegeven, maar kan je daarbij ook ‘(link) functies’ maken. Hieronder een voorbeeld wat je bijv in restriktor kan:
   myConstraints <- ’
   m1 := .Intercept.
   m2 := .Intercept. + facialBurns
   m3 := .Intercept. + gender
   m4 := .Intercept. + facialBurns + gender + gender.facialBurns
   m1 < m4
   m2 < m4
   m3 < m4 ’
   Kan iets soortgelijks ook in de gorica? Of kan het alleen maar obv de namen van de variabelen?

• The function con_weights_boot() of Leonard (so, in restriktor), see: https://rdrr.io/cran/restriktor/man/mixingWeights.html, calculates the weights we also use. We could look into that function and do the same.

Is het mogelijk om het aantal iteraties in het bepalen van de penalty op te geven?
Zo ja, hoe en kan dat in de description/help file worden toegevoegd.
Zo nee, kan je dat realiseren (en dan toevoegen aan help file)?

Rebecca: Ik doe ook iets met ‘combining evidence’. Waarschijnlijk zal ik nog steeds moeten wennen, maar hier lijkt mij het opgeven van restricties in matrices handiger dan in woorden (of anders gezegd: hier lijkt me het opgeven van de parameters voor elk van de studies in de evidence synthesis voorzien van namen in een shiny app erg onhandig). Ik kan dat natuurlijk achter de schermen doen, maar dat is een beetje suf: van restrictiematrix naar hypothese in woorden en dan in jouw deel (dus de gorica R package) weer andersom. Heel misschien moeten we kijken wat we al kunnen inbouwen in de gorica voor wat betreft combining evidence…. Maar dat is weer wat meer (maar komt wel weer een artikel uit met simulaties…).

Laten we ook rel weights zien (ook tov complement)? Ws niet, ms ook een mx maken zoals bain met BF. Dan daar aangeven welke hun max hebben (nl als gelijke LL, dus ook als twee gelijk hebben maat ll niet max is).

Bijgevoegd een file waar je kan zien dat het fout gaat: Ik zal het e-mailen, want ik krijg melding dat deze type file niet gesupport wordt...

Ik krijg zelfs een negatieve penalty voor het complement (en het zou iets dicht bij 12 moeten zijn, iig tussen 0 en 12)...
Ik heb ik restriktor geprobeerd met soortgelijke hypo en dan in een regressie model en daar gaat penalty bepalen voor het complement wel goed; dus ik denk dat het ergens in jouw code zit.... Wellicht moet je code van restriktor aanroepen??

Ook een matrix maken (van relative weights) zoals bain met BFmatrix.
Dat evt alleen opvraagbaar, maar mag ook standaard output zijn.

Evt mbt "Ik wil bijv, als lik = max lik is, aangeven wat max relative weight is of zeggen dat er nu een max is tov hypo’s met zelfde lik.":
Dan in rel weight matrix aangeven welke hun max hebben (nl als gelijke LL, dus ook als twee gelijk hebben maar ll niet max is).