Adding one-sided (calibrated) p-values to results

OHDSI · Nov 13, 2023 · 3560f86 · 3560f86
1 parent f953bad
commit 3560f86
Show file tree

Hide file tree

Showing 5 changed files with 23 additions and 218 deletions.
diff --git a/R/Export.R b/R/Export.R
@@ -515,6 +515,7 @@ exportCohortMethodResults <- function(outputFolder,
       "ci95Lb",
       "ci95Ub",
       "p",
+      "oneSidedP",
       "targetSubjects",
       "comparatorSubjects",
       "targetDays",
@@ -528,6 +529,7 @@ exportCohortMethodResults <- function(outputFolder,
       "calibratedCi95Lb",
       "calibratedCi95Ub",
       "calibratedP",
+      "calibratedOneSidedP",
       "calibratedLogRr",
       "calibratedSeLogRr",
       "targetEstimator"

diff --git a/R/RunAnalyses.R b/R/RunAnalyses.R
@@ -1702,8 +1702,12 @@ summarizeResults <- function(referenceTable, outputFolder, mainFileName, interac
     if (is.null(coefficient)) {
       p <- NA
     } else {
-      z <- coefficient / outcomeModel$outcomeModelTreatmentEstimate$seLogRr
-      p <- 2 * pmin(pnorm(z), 1 - pnorm(z))
+      p <- EmpiricalCalibration::computeTraditionalP(logRr = estimate$logRr,
+                                                     seLogRr = estimate$seLogRr)
+      oneSidedP <- EmpiricalCalibration::computeTraditionalP(logRr = estimate$logRr,
+                                                             seLogRr = estimate$seLogRr,
+                                                             twoSided = FALSE,
+                                                             upper = TRUE)
     }
     pTarget <- outcomeModel$populationCounts$targetExposures /
       (outcomeModel$populationCounts$targetExposures + outcomeModel$populationCounts$comparatorExposures)
@@ -1744,6 +1748,7 @@ summarizeResults <- function(referenceTable, outputFolder, mainFileName, interac
         ci95Lb = if (is.null(coefficient)) NA else exp(ci[1]),
         ci95Ub = if (is.null(coefficient)) NA else exp(ci[2]),
         p = !!p,
+        oneSidedP = !!oneSidedP,
         logRr = if (is.null(coefficient)) NA else coefficient,
         seLogRr = if (is.null(coefficient)) NA else outcomeModel$outcomeModelTreatmentEstimate$seLogRr,
         llr = if (is.null(coefficient)) NA else outcomeModel$outcomeModelTreatmentEstimate$llr,
@@ -1814,7 +1819,14 @@ calibrateGroup <- function(group) {
   if (nrow(ncs) >= 5) {
     null <- EmpiricalCalibration::fitMcmcNull(logRr = ncs$logRr, seLogRr = ncs$seLogRr)
     ease <- EmpiricalCalibration::computeExpectedAbsoluteSystematicError(null)
-    calibratedP <- EmpiricalCalibration::calibrateP(null = null, logRr = group$logRr, seLogRr = group$seLogRr)
+    calibratedP <- EmpiricalCalibration::calibrateP(null = null,
+                                                    logRr = group$logRr,
+                                                    seLogRr = group$seLogRr)
+    calibratedOneSidedP <- EmpiricalCalibration::calibrateP(null = null,
+                                                            logRr = group$logRr,
+                                                            seLogRr = group$seLogRr,
+                                                            twoSided = FALSE,
+                                                            upper = TRUE)
     if (nrow(pcs) >= 5) {
       model <- EmpiricalCalibration::fitSystematicErrorModel(
         logRr = c(ncs$logRr, pcs$logRr),
@@ -1830,6 +1842,7 @@ calibrateGroup <- function(group) {
     group$calibratedCi95Lb <- exp(calibratedCi$logLb95Rr)
     group$calibratedCi95Ub <- exp(calibratedCi$logUb95Rr)
     group$calibratedP <- calibratedP$p
+    group$calibratedOneSidedP <- calibratedOneSidedP$p
     group$calibratedLogRr <- calibratedCi$logRr
     group$calibratedSeLogRr <- calibratedCi$seLogRr
     group$ease <- ease$ease
@@ -1838,6 +1851,7 @@ calibrateGroup <- function(group) {
     group$calibratedCi95Lb <- NA
     group$calibratedCi95Ub <- NA
     group$calibratedP <- NA
+    group$calibratedOneSidedP <- NA
     group$calibratedLogRr <- NA
     group$calibratedSeLogRr <- NA
     group$ease <- NA

diff --git a/inst/csv/resultsDataModelSpecification.csv b/inst/csv/resultsDataModelSpecification.csv
@@ -42,6 +42,7 @@ cm_result,rr,float,No,No,No,No,The estimated relative risk (e.g. the hazard rati
 cm_result,ci_95_lb,float,No,No,No,No,The lower bound of the 95% confidence interval of the relative risk.
 cm_result,ci_95_ub,float,No,No,No,No,The upper bound of the 95% confidence interval of the relative risk.
 cm_result,p,float,No,No,No,No,The two-sided p-value considering the null hypothesis of no effect.
+cm_result,one_sided_p,float,No,No,No,No,he one-sided p-value considering the null hypothesis of RR <= 1.
 cm_result,target_subjects,int,Yes,No,Yes,No,The number of subject in the target cohort.
 cm_result,comparator_subjects,int,Yes,No,Yes,No,The number of subject in the comparator cohort.
 cm_result,target_days,int,Yes,No,No,No,The number of days observed in the target cohort.
@@ -55,6 +56,7 @@ cm_result,calibrated_rr,float,No,No,No,No,The calibrated relative risk.
 cm_result,calibrated_ci_95_lb,float,No,No,No,No,The lower bound of the calibrated 95% confidence interval of the relative risk.
 cm_result,calibrated_ci_95_ub,float,No,No,No,No,The upper bound of the calibrated 95% confidence interval of the relative risk.
 cm_result,calibrated_p,float,No,No,No,No,The calibrated two-sided p-value.
+cm_result,calibrated_one-sided_p,float,No,No,No,No,The calibrated one-sided p-value considering the null hypothesis of RR <= 1.
 cm_result,calibrated_log_rr,float,No,No,No,No,The log of the calibrated relative risk.
 cm_result,calibrated_se_log_rr,float,No,No,No,No,The standard error of the log of the calibrated relative risk.
 cm_result,target_estimator,varchar,Yes,No,No,No,"The target estimator, for example ""att"", ""ate"", ""atu"" or ""ato""."

diff --git a/inst/sql/sql_server/migrations/Migration_1-v5_2_0.sql b/inst/sql/sql_server/migrations/Migration_1-v5_2_0.sql
@@ -1,6 +1,6 @@
 -- Database migrations for version 5.2.0
 -- This migration updates the schema:
--- 1. Table cm_result add target_estimator field
+-- 1. Table cm_result add target_estimator, one_sided_p and calibrated_one_sided_p fields
 -- 2. Tablec m_interaction_result add target_estimator field
 -- 3. Table cm_covariate_balance and cm_shared_covariate_balance add:
 --   - mean_before
@@ -14,8 +14,7 @@
 --    Add:
 --   - generalizability_max_sdm
 --   - generalizability_diagnostic
-
-ALTER TABLE @database_schema.@table_prefixcm_result ADD target_estimator VARCHAR(3);
+ALTER TABLE @database_schema.@table_prefixcm_result ADD target_estimator VARCHAR(3), one_sided_p FLOAT, calibrated_one_sided_p FLOAT;
 ALTER TABLE @database_schema.@table_prefixcm_interaction_result ADD target_estimator VARCHAR(3);
 ALTER TABLE @database_schema.@table_prefixcm_covariate_balance ADD mean_before FLOAT, mean_after FLOAT, target_std_diff FLOAT, comparator_std_diff FLOAT, target_comparator_std_diff FLOAT;
 ALTER TABLE @database_schema.@table_prefixcm_shared_covariate_balance ADD mean_before FLOAT, mean_after FLOAT, target_std_diff FLOAT, comparator_std_diff FLOAT, target_comparator_std_diff FLOAT;

diff --git a/tests/testthat/setup.R b/tests/testthat/setup.R
@@ -1,221 +1,9 @@
 library(Eunomia)
-library(CohortMethod)
-
 connectionDetails <- getEunomiaConnectionDetails()
 Eunomia::createCohorts(connectionDetails)
 
-# fitOutcomeModel ----
-## Study Population ----
-nsaids <- c(1118084, 1124300)
-
-covSettings <- createDefaultCovariateSettings(
-  excludedCovariateConceptIds = nsaids,
-  addDescendantsToExclude = TRUE
-)
-
-sCohortMethodData <- getDbCohortMethodData(
-  connectionDetails = connectionDetails,
-  cdmDatabaseSchema = "main",
-  targetId = 1,
-  comparatorId = 2,
-  outcomeIds = c(3, 4),
-  exposureDatabaseSchema = "main",
-  outcomeDatabaseSchema = "main",
-  exposureTable = "cohort",
-  outcomeTable = "cohort",
-  covariateSettings = covSettings
-)
-
-studyPop <- createStudyPopulation(
-  cohortMethodData = sCohortMethodData,
-  outcomeId = 3,
-  riskWindowEnd = 99999
-)
-
-## CohortMethod Data ----
-sCohortMethodData <- getDbCohortMethodData(
-  connectionDetails = connectionDetails,
-  cdmDatabaseSchema = "main",
-  targetId = 1,
-  comparatorId = 2,
-  outcomeIds = c(3, 4),
-  exposureDatabaseSchema = "main",
-  outcomeDatabaseSchema = "main",
-  exposureTable = "cohort",
-  outcomeTable = "cohort",
-  covariateSettings = covSettings
-)
-
-ps <- createPs(cohortMethodData = sCohortMethodData, population = studyPop)
-studyPopStratisfied <- stratifyByPs(ps, 5)
-studyPopMatched <- matchOnPs(population = ps)
-
-## runCmAnalyses ----
-outputFolder <- tempfile(pattern = "cmData")
-
-covarSettings <- createDefaultCovariateSettings(addDescendantsToExclude = TRUE)
-
-getDbCmDataArgs <- createGetDbCohortMethodDataArgs(
-  washoutPeriod = 183,
-  firstExposureOnly = TRUE,
-  removeDuplicateSubjects = "remove all",
-  covariateSettings = covarSettings
-)
-
-createPsArgs <- createCreatePsArgs(
-  prior = createPrior("laplace", variance = 0.01),
-  estimator = "att"
-)
-
-matchOnPsArgs <- createMatchOnPsArgs(maxRatio = 100)
-
-computeSharedCovBalArgs <- createComputeCovariateBalanceArgs()
-
-computeCovBalArgs <- createComputeCovariateBalanceArgs(
-  covariateFilter = FeatureExtraction::getDefaultTable1Specifications()
-)
-
-truncateIptwArgs <- createTruncateIptwArgs(maxWeight = 10)
-
-tcos1 <- createTargetComparatorOutcomes(
-  targetId = 1,
-  comparatorId = 2,
-  outcomes = list(
-    createOutcome(
-      outcomeId = 3,
-      priorOutcomeLookback = 30
-    ),
-    createOutcome(
-      outcomeId = 4,
-      outcomeOfInterest = FALSE,
-      trueEffectSize = 1
-    )
-  ),
-  excludedCovariateConceptIds = c(1118084, 1124300)
-)
-# Empty cohorts:
-tcos2 <- createTargetComparatorOutcomes(
-  targetId = 998,
-  comparatorId = 999,
-  outcomes = list(
-    createOutcome(
-      outcomeId = 3,
-      priorOutcomeLookback = 30
-    ),
-    createOutcome(
-      outcomeId = 4,
-      outcomeOfInterest = FALSE,
-      trueEffectSize = 1
-    )
-  )
-)
-
-targetComparatorOutcomesList <- list(tcos1, tcos2)
-
-analysesToExclude <- data.frame(
-  targetId = c(998, 998),
-  analysisId = c(3, 4)
-)
-
-createStudyPopArgs1 <- createCreateStudyPopulationArgs(
-  removeSubjectsWithPriorOutcome = TRUE,
-  firstExposureOnly = TRUE,
-  restrictToCommonPeriod = TRUE,
-  removeDuplicateSubjects = "remove all",
-  washoutPeriod = 183,
-  censorAtNewRiskWindow = TRUE,
-  minDaysAtRisk = 1,
-  riskWindowStart = 0,
-  startAnchor = "cohort start",
-  riskWindowEnd = 30,
-  endAnchor = "cohort end"
-)
-
-createStudyPopArgs2 <- createCreateStudyPopulationArgs(
-  removeSubjectsWithPriorOutcome = TRUE,
-  firstExposureOnly = TRUE,
-  restrictToCommonPeriod = TRUE,
-  removeDuplicateSubjects = "keep first",
-  washoutPeriod = 183,
-  censorAtNewRiskWindow = TRUE,
-  minDaysAtRisk = 1,
-  riskWindowStart = 0,
-  startAnchor = "cohort start",
-  riskWindowEnd = 30,
-  endAnchor = "cohort end"
-)
-
-### Analysis 1 ----
-fitOutcomeModelArgs1 <- createFitOutcomeModelArgs(
-  modelType = "cox"
-)
-
-cmAnalysis1 <- createCmAnalysis(
-  analysisId = 1,
-  description = "No matching, simple outcome model",
-  getDbCohortMethodDataArgs = getDbCmDataArgs,
-  createStudyPopArgs = createStudyPopArgs1,
-  fitOutcomeModelArgs = fitOutcomeModelArgs1
-)
-
-### Analysis 2 ----
-fitOutcomeModelArgs2 <- createFitOutcomeModelArgs(
-  modelType = "cox",
-  stratified = TRUE
-)
-
-cmAnalysis2 <- createCmAnalysis(
-  analysisId = 2,
-  description = "Matching",
-  getDbCohortMethodDataArgs = getDbCmDataArgs,
-  createStudyPopArgs = createStudyPopArgs2,
-  createPsArgs = createPsArgs,
-  matchOnPsArgs = matchOnPsArgs,
-  computeSharedCovariateBalanceArgs = computeSharedCovBalArgs,
-  computeCovariateBalanceArgs = computeCovBalArgs,
-  fitOutcomeModelArgs = fitOutcomeModelArgs2
-)
-
-### Analysis 3 ----
-fitOutcomeModelArgs3 <- createFitOutcomeModelArgs(
-  modelType = "cox",
-  inversePtWeighting = TRUE
-)
-cmAnalysis3 <- createCmAnalysis(
-  analysisId = 3,
-  description = "IPTW",
-  getDbCohortMethodDataArgs = getDbCmDataArgs,
-  createStudyPopArgs = createStudyPopArgs2,
-  createPsArgs = createPsArgs,
-  truncateIptwArgs = truncateIptwArgs,
-  computeSharedCovariateBalanceArgs = computeSharedCovBalArgs,
-  fitOutcomeModelArgs = fitOutcomeModelArgs3
-)
-
-### Analysis 4 ----
-fitOutcomeModelArgs4 <- createFitOutcomeModelArgs(
-  modelType = "cox",
-  stratified = TRUE,
-  interactionCovariateIds = 8532001
-)
-
-cmAnalysis4 <- createCmAnalysis(
-  analysisId = 4,
-  description = "Matching with gender interaction",
-  getDbCohortMethodDataArgs = getDbCmDataArgs,
-  createStudyPopArgs = createStudyPopArgs2,
-  createPsArgs = createPsArgs,
-  matchOnPsArgs = matchOnPsArgs,
-  fitOutcomeModelArgs = fitOutcomeModelArgs4
-)
-
-cmAnalysisList <- list(cmAnalysis1, cmAnalysis2, cmAnalysis3, cmAnalysis4)
-
-## Clean-up ----
 withr::defer(
   {
-    unlink(outputFolder)
-
     # Remove the Eunomia database:
     unlink(connectionDetails$server())
     if (getOption("use.devtools.sql_shim", FALSE)) {