R Under development (unstable) (2025-03-11 r87944 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2025 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(qgcomp)
> library(qgcompint)
> set.seed(23)
> dat2 <- simdata_quantized_emm(
+   outcometype="logistic",
+   # sample size
+   n = 100,
+   # correlation between x1 and x2, x3, ...
+   corr=c(0.8, 0.6, 0.3, -0.3, -0.3, -0.3),
+   # model intercept
+   b0=-2,
+   # linear model coefficients for x1, x2, ... at referent level of interacting variable
+   mainterms=c(0.3, -0.1, 0.1, 0.0, 0.3, 0.1, 0.1),
+   # linear model coefficients for product terms between x1, x2, ... and interacting variable
+   prodterms = c(1.0, 0.0, 0.0, 0.0, 0.2, 0.2, 0.2),
+   # type of interacting variable
+   ztype = "categorical",
+   # number of levels of exposure
+   q = 4,
+   # residual variance of y
+   yscale = 2.0
+ )
Warning message:
In .dgm_quantized_logistic_emm(N = n, b0 = b0, mainterms = mainterms,  :
  model implies > 10% of observations with very high/low (<0.5%) outcome probability, which may distort estimates
> 
> head(dat2)
  z x1 x2 x3 x4 x5 x6 x7 y
1 0  1  1  0  3  0  3  2 0
2 2  2  3  2  1  3  1  1 1
3 0  3  0  0  3  0  0  2 1
4 2  1  1  1  0  3  2  2 1
5 2  2  2  2  2  0  0  1 1
6 0  0  0  0  0  1  3  3 0
> table(dat2$z)

 0  1  2 
33 33 34 
> table(dat2$y)

 0  1 
35 65 
> dat2$z = as.factor(dat2$z)
> qfit2 <- qgcomp.emm.glm.noboot(y~x1,
+                            data = dat2,
+                            expnms = paste0("x", 1:1),
+                            emmvar = "z",
+                            q = 4)
> qfit2
## Qgcomp weights/partial effects at z = 0
Scaled effect size (positive direction, sum of positive effects = 0.112)
x1 
 1 
Scaled effect size (negative direction, sum of negative effects = 0)
None


## Mixture slope parameters (delta method CI):

            Estimate Std. Error  Lower CI Upper CI t value  Pr(>|t|)
(Intercept) 0.038714   0.096934 -0.151274  0.22870  0.3994   0.69051
psi1        0.112205   0.049806  0.014587  0.20982  2.2529   0.02657
z1          0.573912   0.138179  0.303086  0.84474  4.1534 7.155e-05
z1:mixture  0.043462   0.072059 -0.097771  0.18470  0.6031   0.54785
z2          0.615288   0.138944  0.342962  0.88761  4.4283 2.539e-05
z2:mixture  0.046243   0.075029 -0.100811  0.19330  0.6163   0.53915
> qfit2$fit

Call:  glm(formula = newform, data = qdata, weights = weights)

Coefficients:
(Intercept)           x1           z1           z2        x1:z1        x1:z2  
    0.03871      0.11220      0.57391      0.61529      0.04346      0.04624  

Degrees of Freedom: 99 Total (i.e. Null);  94 Residual
Null Deviance:	    22.75 
Residual Deviance: 10.77 	AIC: 74.93
> 
> proc.time()
   user  system elapsed 
   2.56    0.53    3.07