Learning individualized treatment regimes (ITR) using observational data holds great interest in various fields, as treatment recommendations based on individual characteristics may improve individual treatment benefits with a reduced cost. It has long been observed that different individuals may respond to a certain treatment with significant heterogeneity. ITR can be defined as a mapping between individual characteristics to a treatment assignment. The optimal ITR is the treatment assignment that maximizes expected individual treatment effects. Rooted from personalized medicine, many studies and applications of ITR are in medical fields and clinical practice. Heterogeneous responses are also well documented in educational interventions. However, unlike the efficacy study in medical studies, educational interventions are often not randomized. Study results often suffer greatly from self-selection bias. Besides the intervention itself, the efficacy and effectiveness of interventions usually interact with a wide range of confounders. In this study, we propose a novel algorithm to extend random forest of interaction trees to Casual Effect Random Forest of Interaction Trees (CERFIT) for learning individualized treatment effects and regimes. We first consider the study under a binary treatment setting. Each interaction tree recursively partitions the data into two subgroups with greatest heterogeneity of treatment effect. By integrating propensity score into the tree growing process, subgroups from the proposed CERFIT not only have maximized treatment effect differences, but also similar baseline covariates. Thus it allows for the estimation of the individualized treatment effects using observational data. In addition, we also propose to use residuals from linear models instead of the original responses in the algorithm. By doing so, the numerical stability of the algorithm is greatly improved, which leads to an improved prediction accuracy. We then consider the learning problem under non-binary treatment settings. For multiple treatments, through recursively partitioning data into two subgroups with greatest treatment effects heterogeneity with respect to two randomly selected treatment groups, the algorithm transforms the multiple learning ITR into a binary task. Similarly, continuous treatment can be handled through recursively partitioning the data into subgroups with greatest homogeneity in terms of the association between the response and the treatment within a child node. For all treatment settings, the CERFIT provides variable importance ranking in terms of treatment effects. Extensive simulation studies for assessing estimation accuracy and variable importance ranking are presented. CERFIT demonstrates competitive performance among all competing methods in simulation studies. The methods are also illustrated through an assessment of a voluntary education intervention for binary treatment setting and learning optimal ITR among multiple interventions for non-binary treatments using data from a large public university.