R Under development (unstable) (2024-10-11 r87227 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2024 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.

> #rm(list=ls())
> # save(edurose_mediation_20181126,file = "data/edurose_mediation_20181126.rda",compress = "xz")
> 
> #library("htetree")
> # load("data/edurose_mediation_20181126.rda")
> 
> # construct the simulated data based on Athey's data
> #install.packages("data.table",repos = NULL,
> #                 type = "source")
> #library("data.table")
> #install.packages("causalTree",
> #                 repos = "https://jiahui1902.github.io/drat/",
> #                 type = "source")
> # install.packages("htetree",
> #                  repos = "https://jiahui1902.github.io/drat/",
> #                  type = "source")
> # library(causalTree)
> 
> library(rpart)
> library(htetree)
> hte_causalTree(outcomevariable="outcome",
+                data=data.frame("confounder"=c(0, 1, 1, 0, 1, 1),
+                                "treatment"=c(0,0,0,1,1,1),
+                                "prop_score"=c(0.4, 0.4, 0.5, 0.6, 0.6, 0.7),
+                                "outcome"=c(1, 2, 2, 1, 4, 4)),
+                treatment_indicator = "treatment",
+                ps_indicator = "prop_score", covariates = "confounder")
$predictedTE
  outcome_predictedTE
1            1.333333
2            1.333333
3            1.333333
4            1.333333
5            1.333333
6            1.333333

$tree
n= 6 

node), split, n, deviance, yval
      * denotes terminal node

1) root 6 383.1111 1.333333 *

> # causalTree(y~ x1 + x2 + x3 + x4, data = simulation.1,
> #     treatment = simulation.1$treatment,
> #     split.Rule = "CT", cv.option = "CT", split.Honest = TRUE, cv.Honest = TRUE,
> #     split.Bucket = F, xval = 5,
> #     cp = 0, minsize = 20, propensity = 0.5)
> data("simulation.1")
> # estimate the propensity score
> fit <- glm(treatment~x1+x2+x3+x4+x5+x6+x7+x8+x9+x10,
+            data=simulation.1,
+            family = "binomial")
> simulation.1$ps_score <- predict(fit,type = "response")
> linear_terms <- paste0("x",1:10)
> 
> # estimate the model with our package
> set.seed(1)
> lb <- c(paste0("var",1:10),"propensity score")
> names(lb) <- c(paste0("x",1:10),"ps_score")
> 
> fit_drawplot <- htetree::hte_ipw(outcomevariable = 'y',
+                minsize=20,crossvalidation = 40,negative = TRUE,
+                data = simulation.1,
+                ps_indicator = "ps_score",
+                covariates = c(linear_terms, "ps_score"),
+                drawplot = TRUE,treatment_indicator = "treatment",
+                # no_indicater = '_IPW_simulation',
+                legend.x = 0.1,legend.y = 0.25,varlabel = lb)
box.palette (not diverging): #0000FF (blue) to #DAA520 (goldenrod)
cex 1   xlim c(0, 1)   ylim c(0, 1)
> 
> fit_noplot <- htetree::hte_ipw(outcomevariable = 'y',
+                minsize=20,crossvalidation = 40,negative = TRUE,
+                data = simulation.1,
+                ps_indicator = "ps_score",
+                covariates = c(linear_terms, "ps_score"),
+                drawplot = TRUE,treatment_indicator = "treatment",
+                # no_indicater = '_IPW_simulation',
+                legend.x = 0.1,legend.y = 0.25,varlabel = lb)
box.palette (not diverging): #0000FF (blue) to #DAA520 (goldenrod)
cex 1   xlim c(0, 1)   ylim c(0, 1)
> 
> # hte_plot_line(model = xxx,data = simulation.1,
> #               treatment_indicator = "treatment",
> #               outcomevariable = 'y',
> #               propensity_score = "ps_score",gamma = 0.5,lambda = 0.5)
> # hte_plot_line(model = xxx,data = simulation.1,
> #               treatment_indicator = "treatment",
> #               outcomevariable = 'y',
> #               propensity_score = "ps_score")
> 
> # library(htetree)
> # hte_plot(model = fit_drawplot,data = simulation.1,
> #          treatment_indicator = "treatment",
> #          outcomevariable = 'y',
> #          propensity_score = "ps_score")
> #
> # hte_plot_line(model = fit_drawplot,data = simulation.1,
> #               treatment_indicator = "treatment",
> #               outcomevariable = 'y',
> #               propensity_score = "ps_score")
> 
> 
> # ps_indicator = 'ps_score'
> # covs <- c("x1", "x2", "x3",
> #           "x4", "x5", "x6", "x7",
> #           "x8", "x9", "x10", "ps_score")
> # xxx2 <- hte_matchinginleaves(outcomevariable = 'y',
> #                              data = simulation.1,
> #                              drawplot = TRUE,
> #                              ps_indicator = "ps_score",
> #                              treatment_indicator = "treatment",
> #                              covariates=covs,
> #                              con.num=4)
> 
> ## Dynamic visualization for IPW model (using Shiny)
> # THIS IS NOT SHOWING UP (shows a blank page)
> # The runDynamic function runs the visualization without saving any of the files
> # runDynamic(fit_drawplot, simulation.1,
> #            outcomevariable = "y", treatment_indicator = "treatment",
> #            propensity_score="ps_score")
> 
> # The files for runDynamic are saved in the temporary directory
> # The files can be cleared manually using the clearTemp() function, or will automatically be cleared when you close R
> # clearTemp()
> 
> proc.time()
   user  system elapsed 
   2.04    0.32    2.34