fix anova type 3

[clarkliming]

close #459

[clarkliming]

after #455 merged we may further improve the code

[danielinteractive]

thanks @clarkliming , I also merged the other PR, so can go from updated version now

[danielinteractive]

Thanks @clarkliming

[github-actions]

Code Coverage Summary

Filename                    Stmts    Miss  Cover    Missing
------------------------  -------  ------  -------  ----------------------------
R/between-within.R             59       0  100.00%
R/component.R                  69       0  100.00%
R/cov_struct.R                 97       1  98.97%   407
R/empirical.R                   7       0  100.00%
R/fit.R                       229       3  98.69%   420, 481, 511
R/interop-car.R               130       1  99.23%   9
R/interop-emmeans.R            51       0  100.00%
R/interop-parsnip.R            59       1  98.31%   12
R/kenwardroger.R               92       2  97.83%   41, 63
R/mmrm-methods.R              140       0  100.00%
R/residual.R                    8       0  100.00%
R/satterthwaite.R             116      12  89.66%   238-249
R/skipping.R                    8       0  100.00%
R/testing.R                    64       4  93.75%   29, 31, 80, 82
R/tidiers.R                    72       2  97.22%   46-47
R/tmb-methods.R               283       3  98.94%   279-280, 336
R/tmb.R                       295       0  100.00%
R/utils-formula.R              27       0  100.00%
R/utils-nse.R                  16       0  100.00%
R/utils.R                     199      12  93.97%   279-289, 437, 466
R/zzz.R                        70      24  65.71%   7-22, 55-60, 90, 118, 138
src/chol_cache.h               63       0  100.00%
src/covariance.h              101       1  99.01%   177
src/derivatives.h             126       0  100.00%
src/empirical.cpp              72       0  100.00%
src/exports.cpp                47       0  100.00%
src/jacobian.cpp               47       1  97.87%   54
src/kr_comp.cpp                56       0  100.00%
src/mmrm.cpp                   76       0  100.00%
src/predict.cpp                93       0  100.00%
src/test-chol_cache.cpp        58       5  91.38%   9, 18, 26, 55, 62
src/test-covariance.cpp       123       5  95.93%   9, 29, 40, 61, 72
src/test-derivatives.cpp      108       7  93.52%   44, 53, 62, 85, 94, 106, 124
src/test-utils.cpp            195       7  96.41%   9, 16, 24, 34, 44, 57, 119
src/testthat-helpers.h         15       5  66.67%   36-37, 41, 50, 53
src/utils.h                    78       0  100.00%
TOTAL                        3349      96  97.13%

Diff against main

Filename           Stmts    Miss  Cover
---------------  -------  ------  -------
R/interop-car.R      +58      -2  +3.40%
TOTAL                +58      -2  +0.11%

Results for commit: 1bda930

Minimum allowed coverage is 80%

♻️ This comment has been updated with latest results

[github-actions]

Unit Tests Summary

    1 files     46 suites   27s ⏱️
  505 tests   465 ✅ 40 💤 0 ❌
1 954 runs  1 908 ✅ 46 💤 0 ❌

Results for commit 1bda930.

♻️ This comment has been updated with latest results.

[github-actions]

Unit Test Performance Difference

Additional test case details

Test Suite	$Status$	Time on main	$±Time$	Test Case
car	👶		$+0.01$	h_first_contain_categorical_works_as_expected
car	👶		$+0.26$	h_get_contrast_works_even_if_only_interaction_term_exists
car	👶		$+0.18$	h_get_contrast_works_even_if_the_interaction_term_order_changes
car	👶		$+0.18$	h_get_contrast_works_for_higher_order_interaction
car	👶		$+0.06$	h_get_contrast_works_if_intercept_is_not_given
car	👶		$+0.01$	h_get_index_works_as_expected
car	👶		$+0.01$	h_obtain_lvls_works_as_expected

Results for commit 61873e5

♻️ This comment has been updated with latest results.

[danielinteractive]

Thanks @clarkliming , just a few minor suggestions

[InkaRo]

Thank you for the quick fix! At first glance, it works fine now. I wanted to run some more thorough checks and also have a look at 3-way interactions. Thanks again!

[kkmann]

It seems that for the 3 way case, there is still a discrepancy.

> library(mmrm)
> library(car)
> 
> mod <- mmrm(
+     formula = FEV1 ~ ARMCD + RACE + AVISIT + RACE * AVISIT * ARMCD + FEV1_BL + us(AVISIT | USUBJID),
+     data = fev_data,
+     control = mmrm_control(
+         method = "Kenward-Roger",
+         vcov = "Kenward-Roger-Linear"
+       ),
+     reml = TRUE
+   )
> 
> car::Anova(mod, type = "III") |>
+   print()
Analysis of Fixed Effect Table (Type III F tests)
                  Num Df Denom Df F Statistic   Pr(>=F)    
ARMCD                  1   169.24      32.726 4.712e-08 ***
RACE                   2   168.43      24.655 4.058e-10 ***
AVISIT                 3   155.86     129.283 < 2.2e-16 ***
FEV1_BL                1   187.76      31.571 6.856e-08 ***
RACE:AVISIT            6   213.76       1.754    0.1099    
ARMCD:RACE             2   137.39       0.199    0.8201    
ARMCD:AVISIT           3   156.17       1.685    0.1725    
ARMCD:RACE:AVISIT      6   198.34       1.269    0.2732    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

SAS reports

The difference is in the two-way interactions now, 3-way and marginals are fine.

I must admit that I am struggling a bit with the code in h_get_contrast() without any comments it is quite dense. Maybe it would help to flesh out the contrast definition in the documentation first?

[clarkliming]

Hi @kkmann I updated again, hopefully this cover all cases

[kkmann]

@InkaRo could you share the current discrepancy case directly in this PR so that we can discuss it here?

@clarkliming it seems that the contrast matrices are stil la bit off in the three way case. I am sorry, but it is quite difficult for Inka and me to reverse engineer the h_* functions without more specs / comments.

Could we maybe document hiw we expect the the h_get_contrast() function to behave and test it directly for higher order interactions?

mmrm/R/interop-car.R

Lines 17 to 28 in 875322e

	#' Obtain Contrast for Specified Effect
	#'
	#' This is support function to obtain contrast matrix for type II/III testing.
	#'
	#' @param object (`mmrm`)\cr the fitted MMRM.
	#' @param effect (`string`) the name of the effect.
	#' @param type (`string`) type of test, "II", "III", '2', or '3'.
	#' @param tol (`numeric`) threshold blow which values are treated as 0.
	#'
	#' @return A `matrix` of the contrast.
	#'
	#' @keywords internal

[clarkliming]

hi @kkmann the details can be found at https://openpharma.github.io/mmrm/latest-tag/articles/hypothesis_testing.html; and this vignettes is just following the SAS. So in what cases do you observe a difference in the contrast?

[danielinteractive]

@kkmann @InkaRo @clarkliming From the sidelines: Maybe it would be easier to have a short meeting to clarify on this topic?

[InkaRo]

A meeting is a good idea. Do you want to set sth up @danielinteractive?

[danielinteractive]

@InkaRo we actually have one today at 2.30 pm CET if you have time? You could let me know your email via [email protected] and I can add you

[kkmann]

Hi @clarkliming, @InkaRo and I just double checked on the latest published commit 875322e and the problem persists. The problem seems to be that the contrasts still need to be sorted in the order of appearance in the formula statement.

[InkaRo]

@clarkliming @danielinteractive sorry, my bad, I did not use load_all(), so I still had the old functions. It works fine now, contrasts look good.

[kkmann]

Apologies from my end as well. Missed it to. Thanks @clarkliming for seeing this through - hope we get that fix to CRAN soon x)

[kkmann]

🚀 works for us, thanks so much for addressing thsi so quickly @clarkliming

There seems to be an issue with a fa icon on the website, remove?