Error when running doIMRMC on Viper data set

[marshhar]

Hi, I'm interested in using the iMRMC R package to analyze data for a non-inferiority trial I'm running. As a first step, I have been trying to replicate the analysis here: https://github.com/DIDSR/iMRMC/wiki/Non-inferiority-by-iMRMC.

However, when I execute this line of code:

results <- doIMRMC(dfMRMC_AI)

I get this error:
Error in if (covAUC < 0) { : missing value where TRUE/FALSE needed

I also tried my own dataset and got the same exception. Do you have any thoughts on what I could be doing wrong? Apologies if this answer is obvious, but I'm new to R.

[marshhar]

I've done a bit of digging and this is the best I've come up with so far.

The error occurs in the doAUCmrmcCov function. Specifically, in the lines of code below, all values in both numer1234 and denom1234 are set to 0:

numer1234 <- colSums(perReaderPair.1[index.TF, c("numer1", "numer2", "numer3", "numer4")])

denom1234 <- colSums(perReaderPair.1[index.TF, c("denom1", "denom2", "denom3", "denom4")])

Later, the former is divided by the latter:

m.biased <- numer/denom

Which results in NaN for the first four values of m.biased. This eventually propagates to covAUC, also set to NaN.

covAUC <- sum(m.coeff * m)

This is later used in a conditional statement:

if (covAUC < 0)

Which throws the exception.

I'm afraid I am not familiar enough with this field to know what is supposed to be occurring. Should the first 4 entries of m.biased be 0 instead of NaN? Any help would be greatly appreciated!

[brandon-gallas]

I believe you have a problem formatting your data. Without the data, I am unable to diagnose your problem. I am also just checking in today and will be on leave. Additional responses will be delayed.

[marshhar]

Thanks for your response, and for the heads up re: delayed responses going forward.

Regarding data, it came from installing this package (viperData 1.0.2): https://github.com/DIDSR/viperData/releases/download/1.0.2/viperData_1.0.2.tar.gz

I'm using iMRMC 2.1.0

The code was copied and pasted from https://github.com/DIDSR/iMRMC/wiki/Non-inferiority-by-iMRMC. For convenience, I've included it all below:

#Load iMRMC package
library(iMRMC)

#Load Viper data 
library(viperData)
viperObs <- viperData::viperObs

#Choose one sub-study => screeningHighP
screeningHighP <- viperObs[viperObs$subStudyLabel == "screeningHighP", , drop = FALSE]
screeningHighP$caseID <- as.character(screeningHighP$caseID)

#Create iMRMC dataframe
dfMRMC <- createIMRMCdf(screeningHighP, 
                        keyColumns = list(readerID = "readerID",
                                          caseID = "caseID",
                                          modalityID = "modalityID",
                                          score = "score",
                                          truth = "cancerStatus455"),
                        truePositiveFactor = "1")

#Subset FFDM reads from all readers
dfMRMC_onlyFFDM <- subset(dfMRMC, modalityID == "FFDM")

#Subset SFM reads from all readers
dfMRMC_onlySFM <- subset(dfMRMC, modalityID == "SFM")

#Subset truth reads from all readers
dfMRMC_onlyTruth <- subset(dfMRMC, modalityID == "truth")

#Create a dataframe for AI simulated reads from FFDM modality
dfMRMC_AImod <- data.frame()

#Loop through all cases in caseID and isolate each case from the chosen readers just once
for(i in unique(dfMRMC_onlyFFDM$readerID)) { 
  remainingCases <- dfMRMC_onlyFFDM[!dfMRMC_onlyFFDM$caseID %in% 
                                      dfMRMC_AImod$caseID, , drop = FALSE]
  dfMRMC_AImod <- rbind(dfMRMC_AImod, subset(remainingCases, readerID == i))
}

#Change modalityID to "AI" for all cases
dfMRMC_AImod$modalityID <- "AI"

#Change readerID to "AI.1" for all cases
dfMRMC_AImod$readerID <- "AI.1"

#Duplicate those cases and rename with readerID "AI.2" 
dfMRMC_AImod2 <- dfMRMC_AImod
dfMRMC_AImod2$readerID <- "AI.2"

#Bind dataframe with 2 AI "readers" with duplicate scores on all cases
dfMRMC_AImod <- rbind(dfMRMC_AImod, dfMRMC_AImod2)

#Bind full dataframe with only SFM, truth and AI rows 
dfMRMC_AI <- rbind(dfMRMC_onlyTruth, dfMRMC_onlySFM, dfMRMC_AImod)

#Run iMRMC analysis function
results <- doIMRMC(dfMRMC_AI)

#Isolate data frame of reader-averaged U-stats based performance results
uStat <- results[['Ustat']] 

#Capture degrees of freedom
#The third row we are calling here assesses the accuracy of modA against modB
dfBDG <- uStat[3, 'dfBDG']

#Capture variance of the difference between modalities, ie. the total variance
varAUCAminusAUCB <- uStat[3, 'varAUCAminusAUCB']

#Define all variables used in the t-statistic calculation 
marginParam = 0.05
aucA = uStat[3, 'AUCA'] 

[marshhar]

Hi, I just wanted to follow up to see if anyone had a chance to look into this.

[brandon-gallas]

I have replicated your error. Thank you for reporting it.

We did test the code in the non-inferiority markdown file. That was last edited in May 2023.

We did a major update that was released Sep. 3, 2024. It seems that the testing did not capture the issue you have found: dividing by 0. This comes from a study design where readers in modality A are not in modality B.

This is a serious bug. I am grateful you logged this issue. I am on it and will keep you updated.

[brandon-gallas]

Hello @marshhar.

I fixed the bug and created two new tests. I also moved the R markdown file into the iMRMC repository for easier access. You can find the update in this branch: https://github.com/DIDSR/iMRMC/tree/debugUnpairedReaders

• New test for unpaired readers: https://github.com/DIDSR/iMRMC/blob/debugUnpairedReaders/Rpackage/iMRMC/tests/testthat/test_doIMRMC_R4.R
• New test for R markdown: https://github.com/DIDSR/iMRMC/blob/debugUnpairedReaders/Rpackage/iMRMC/tests/testthat/test_markdown.R
• R markdown in iMRMC repo: https://github.com/DIDSR/iMRMC/blob/debugUnpairedReaders/Rpackage/iMRMC/inst/extra/000-non-inferiority-by-iMRMC/Non-inferiority-by-iMRMC.rmd

I am asking a colleague to proofread the wiki page. I would be happy if you could check things out too. Thanks, in advance.

Also, I would like to ask you to update your profile to show your real name and affiliation. I like to know who I am helping. Maybe we will meet in real life someday.

All the best. Brandon

[marshhar]

Thank you so much Brandon! I got the Viper example to work.

However, when I try to repurpose the code for my own preliminary data (either at the patient level or lesion leve, both attached), like so:

data <- read.csv(file_path, colClasses = c("character", "character", "character", "numeric", "character"))

dfMRMC <- createIMRMCdf(data, 
                        keyColumns = list(readerID = "readerID",
                                          caseID = "caseID",
                                          modalityID = "modalityID",
                                          score = "score",
                                          truth = "truth"),
                        truePositiveFactor = "1")

results <- doIMRMC(dfMRMC)

I get the following error:

Error in Ustat.1$dfD(Ustat.2$BCK.D.b/Ustat.2$dfD) : 
  attempt to apply non-function
3.
doDFdifference(studySize, Ustat.pair, Ustat.full) at doAUCMRMC.R#397
2.
doAUCmrmc(data, flagROC = TRUE) at iMRMC_R.R#94
1.
doIMRMC(dfMRMC)

It seems to me that the issue starts at line 383 of doAUCMRMC.R:

perReaderPair.12 <- doAUCperReader(mrmcDF.modality.1, mrmcDF.modality.2)

Which returns a dataframe with a value of 0 for numer1, numer2, numer3, denom1, denom2, denom3, and NaN for covAUC, M1, M2, and M3.

Is this similar to the original problem with the Viper dataset? Or is this a problem with my data? (Each subject is randomized to one modality or the other and interpreted by 3 readers; no subject is imaged with both modalities.)

Harry

iMRMC input patient level for R.csv

iMRMC input lesion level for R.csv

[brandon-gallas]

I have an update.

Thanks for sharing data and code. I found a bug. Code was a missing comma related to managing the analysis of your study design; no subject is imaged with both modalities. When I originally fixed your issue where it started (non-inferiority wiki page), the bug was also related to the study design; no reader participated in both modalities. At that time, I thought I fixed the management related to your study design, except for the comma. I had planned on testing this but wanted to respond to you since your issue seemed fixed. Such is life. Regrets.

I took a peek at your results. The patientLevel data is a very small study. Some variance component estimates are negative, so be careful (avoid) working with those. The lesionLevel data has bigger problems. The problem starts with the number of readers. That's not really enough to do this MRMC analysis. I will continue with my tracking and testing to see if there are any other bugs, but you may need to show the per-reader results and not the MRMC results.

Thanks for providing the code and data. That helped a lot. I'll update you again in a day or two.

Brandon

[marshhar]

That's very helpful, thanks for your expertise on this. I did have some questions:

1. I can increase the number of readers - do you have any thoughts on what number might be appropriate to target?
2. I'm surprised the lesion data was more problematic than the patient data, as there were many more observations in the former. Is the issue related to the fact that some lesions were only identified by one or two readers rather than all three?
3. Re: sample size for patient data, the target is n = 200. Is there a way I can check whether this is sufficiently powered?

[brandon-gallas]

My ability to follow through on issues and comments is severely impacted by the recent firings of government issues. Please be patient. Here's some quick responses to tide you over:

1. I can increase the number of readers - do you have any thoughts on what number might be appropriate to target?

   a. 6-8 could be enough. more is always better (says the statistician)

2. I'm surprised the lesion data was more problematic than the patient data, as there were many more observations in the former. Is the issue related to the fact that some lesions were only identified by one or two readers rather than all three?

   a. Yes. That is a likely reason, but I was going to dig in and see.

3. Re: sample size for patient data, the target is n = 200. Is there a way I can check whether this is sufficiently powered?

   a. Yes. We hope to have more tools to help. For now, you can look at this document: https://raw.githubusercontent.com/DIDSR/viperData/master/inst/docs/viperSupplementaryMaterials/viperSupplementaryMaterials.Rmd

[marshhar]

Thanks so much, this is very helpful. Sorry about the issues created by the government firings, hope things improve; I'm in no hurry - looking forward to hearing back when you get a chance.

[marshhar]

Hi Brandon,

I hope things have settled down on your end. I just wanted to send a quick update.

I modified the lesion data to include only lesions identified by all three readers (attached). This adjustment seems to have helped with the issues we discussed — I’ve attached the updated results for both patient-level and lesion-level analyses.

One thing is still puzzling me. To roughly estimate the impact of doubling the number of readers, I duplicated the existing three-reader data, renamed the readers, and added these simulated readers to the dataset (also attached). I expected this to improve the results, but instead, the t-statistic and p-value actually decreased (from 1.058 and 0.294 to 0.916 and 0.361, respectively).

Do you have any insights into why this might be happening?

iMRMC input lesion level for R 6 readers.csv

iMRMC input lesion level for R.csv

Interim Analysis.pdf

[brandon-gallas]

Thanks. This is in my queue, but I don't have time to address. Everything is day-to-day on this kind of work. Regrets.

[marshhar]

I completely understand, I appreciate the help you've provided thus far!