Predictive computational modeling to define effective treatment strategies for bone metastatic prostate cancer

The ability to rapidly assess the efficacy of therapeutic strategies for incurable bone metastatic prostate cancer is an urgent need. Pre-clinical in vivo models are limited in their ability to define the temporal effects of therapies on simultaneous multicellular interactions in the cancer-bone microenvironment. Integrating biological and computational modeling approaches can overcome this limitation. Here, we generated a biologically driven discrete hybrid cellular automaton (HCA) model of bone metastatic prostate cancer to identify the optimal therapeutic window for putative targeted therapies. As proof of principle, we focused on TGFβ because of its known pleiotropic cellular effects. HCA simulations predict an optimal effect for TGFβ inhibition in a pre-metastatic setting with quantitative outputs indicating a significant impact on prostate cancer cell viability, osteoclast formation and osteoblast differentiation. In silico predictions were validated in vivo with models of bone metastatic prostate cancer (PAIII and C4-2B). Analysis of human bone metastatic prostate cancer specimens reveals heterogeneous cancer cell use of TGFβ. Patient specific information was seeded into the HCA model to predict the effect of TGFβ inhibitor treatment on disease evolution. Collectively, we demonstrate how an integrated computational/biological approach can rapidly optimize the efficacy of potential targeted therapies on bone metastatic prostate cancer.

The ability to rapidly assess the efficacy of therapeutic strategies for incurable bone metastatic prostate cancer is an urgent need. Pre-clinical in vivo models are limited in their ability to define the temporal effects of therapies on simultaneous multicellular interactions in the cancer-bone microenvironment. Integrating biological and computational modeling approaches can overcome this limitation. Here, we generated a biologically driven discrete hybrid cellular automaton (HCA) model of bone metastatic prostate cancer to identify the optimal therapeutic window for putative targeted therapies. As proof of principle, we focused on TGFβ because of its known pleiotropic cellular effects. HCA simulations predict an optimal effect for TGFβ inhibition in a pre-metastatic setting with quantitative outputs indicating a significant impact on prostate cancer cell viability, osteoclast formation and osteoblast differentiation. In silico predictions were validated in vivo with models of bone metastatic prostate cancer (PAIII and C4-2B). Analysis of human bone metastatic prostate cancer specimens reveals heterogeneous cancer cell use of TGFβ. Patient specific information was seeded into the HCA model to predict the effect of TGFβ inhibitor treatment on disease evolution. Collectively, we demonstrate how an integrated computational/biological approach can rapidly optimize the efficacy of potential targeted therapies on bone metastatic prostate cancer.
Metastatic castrate resistant prostate cancer (mCRPC) typically manifests in the skeleton and is currently incurable 1,2 . In the bone microenvironment, prostate cancer cells hijack the normal bone remodeling process to create a "vicious cycle" of extensive bone formation and destruction 3 . Key mechanisms facilitating the cross-talk between the cancer and host compartment include the induction of receptor activator of nuclear κ B ligand (RANKL) expression and the release of sequestered growth factors from the bone matrix. Bone is a rich source of transforming growth factorβ (TGFβ ) and the role for this pleiotropic factor in promoting the survival and growth of bone metastatic cancers has been well described 4,5 . The molecular complexity of the circuitry driving this cycle has expanded tremendously in the past two decades revealing many potential targets for therapeutic intervention. The question remains however as to how to translate these potential therapies to the clinic. Biological experimentation and pre-clinical mouse models can be used to define the impact of putative therapies but are limited in their ability to dissect the potential dynamic and simultaneous effects on the multi-cellular tumor-bone microenvironment. One potential alternative approach is the integration of experimentally measured biological parameters with computational models to tackle the multi-scale nature of the disease 6 . Numerous computational models successfully demonstrate the feasibility of the approach [7][8][9][10][11][12][13][14] . Starting from existing experimental or clinical data it is possible to use statistical frameworks such as Approximate Bayesian Computation (ABC) to identify, in a "top-down" manner, the importance of unknown parameters in disease progression by applying a distribution of probability on those factors 15 . Conversely, agent based models, such as discrete-continuum Hybrid Cellular

Results
Computational modeling of TGFβ inhibition in normal bone remodeling and in bone metastatic prostate cancer. TGFβ is known to have concentration dependent pleiotropic effects on osteoblasts and osteoclasts [26][27][28] . In silico, the ability of stromal cells to respond to varying TGFβ concentrations (0.1 to 10 ng/ml) was integrated into our HCA of normal bone remodeling, the bone multicellular unit (BMU) (see Equations 1, 2, and Table S1) 18 . To obtain a realistic readout of the level of TGFβ inhibition that could be achievable in vivo, we treated mice with a TGFβ neutralizing antibody (1D11) at a dose previously used in the literature (10 mg/kg) and consistent with clinical trials performed for a humanized version of the 1D11 TGFβ neutralizing antibody, fresolimumab/GC1008 22,29,30 . We observed that the TGFβ neutralizing antibody significantly reduced circulating TGFβ serum levels by up to 80% compared to IgG control treated animals ( Supplementary Fig. S1). Using phospho-SMAD2 as a surrogate for TGFβ activity 31 , we also observed that the TGFβ neutralizing antibody could inhibit TGFβ activity 50-80% in tumor naïve and tumor bearing tissues (Supplementary Fig. S1 and Fig. 3f). Based on this in vivo information, we applied the TGFβ inhibitor to the normal BMU at a level of 80% in silico. The stochastic nature of the BMU allows for variation and statistical analysis of simulation outputs. The results of multiple simulations (n = 29/group) show that TGFβ inhibition significantly promoted bone formation (9% increase) over a 75-day period by enhancing osteoblast expansion and differentiation while limiting osteoclast viability (Fig. 1a, Supplementary Fig. S1, and Table S2 ). Importantly, these in silico results are consistent with previous in vivo studies and support the robustness of the parameters used to power the computational model 29,32 .
The vicious cycle paradigm suggests that metastatic prostate cancer cells utilize TGFβ signaling to promote their survival and growth. We therefore seeded the computational model with TGFβ responsive prostate cancer cells. Once the vicious cycle was established at day 80, we initiated TGFβ inhibition (post-treatment scenario). Simulations (n = 24/group) revealed that TGFβ inhibition reduced cancer growth by approximately 15%, but only when the inhibitor was applied at a constant 99% level of efficacy until day 250 ( Supplementary Fig. S2 and Video 1). At a more biologically relevant level of 80% inhibition, we observed little difference in cancer cell growth between the control and treatment groups. Surprisingly, and in contrast to the observed effects of TGFβ inhibition on the normal BMU, we also observed no difference in osteogenesis between the control and treatment groups even at later stages ( Supplementary Fig. S2, Table S3, and Video 1). Taken together, these results suggest that the treatment of established and actively growing bone metastases with TGFβ inhibitors would have no impact on the progression of the lesions unless the inhibitor was > 99% effective.
A major advantage of the computational model is that it can be used to explore therapeutic windows of efficacy for putative inhibitors. Simulations (n = 24) revealed that applying the TGFβ inhibitor in silico at day 1 prior to the seeding of the cancer cells (pre-treatment scenario) even at a level of 20% efficacy significantly reduced tumor burden over time by ≥ 65% (Fig 1b, Supplementary Fig. S2, Table S3, and Video 2). Interestingly, TGFβ inhibition resulted in a small but significant increase in cancer-induced osteogenesis compared to control during early tumor progression (Day 100). However, at the end of the simulations there was significantly less cancer-induced osteogenesis in the TGFβ treated group compared to control (Fig. 1b, Supplementary Fig. S2, and These in silico results predict that applying TGFβ inhibitors in a preventative manner will reduce the growth of metastatic prostate cancer without exacerbating cancer induced osteogenesis. In vivo validation of computational model predictions. Analysis of human specimens of bone metastatic prostate cancer derived from deidentified cancer patients at the Moffitt Cancer Center (n = 20) show that TGFβ ligand and receptors are expressed and pSMAD2 staining indicates that the TGFβ pathway is active (Fig. 2a). In vitro analysis of prostate cancer cell lines that can grow in the bone microenvironment identified the PAIII cell line as reflecting the TGFβ receptor and growth factor producing (TRP) status observed in human specimens (Fig. 2b,c). We also noted that the PAIII cell line was sensitive to inhibition with TGFβ inhibitors (Fig. 2d,e). We therefore initially chose this cell line to test computational model predictions.
Dissecting cancer cell behavior in silico illustrates TGFβ inhibition directly limits growth over time by impacting cancer cell proliferation ( Fig. 3a-d). To determine the validity of the computational outputs, 6-week-old male SCID Beige mice were pre-treated with either a TGFβ inhibitor (TGFβi -1D11, 10 mg/Kg, 3× weekly; n = 10) or an isotype control IgG (Control-13C4, 10 mg/Kg 3× weekly; n = 8) and subsequently inoculated with luciferase-expressing PAIII cells. Bioluminescence analysis revealed a significant reduction in tumor growth in the TGFβ inhibitor treated group compared to controls and, as expected, reduced pSMAD2 and AKT phosphorylation (Fig. 3e,f). We further found significant reductions in proliferation (40%) and increases in apoptosis (70%) between the TGFβ inhibitor treated and control groups (Fig. 3g). The HCA model is based on humanized In silico control and TGFβ inhibitor treated simulations in normal (BMU, n = 29/group, a) and bone metastatic prostate cancer (PCA-Bone Mets, n = 24/group, b) scenarios. Representative images of simulation runs at indicated time points are shown with magnified insets. TGFβ inhibitor was applied at day 1 for all simulations (pre-treatment scenario). Cell populations analyzed include mesenchymal stem cell (MSC), osteoblast precursors (pOB), osteoblasts, osteoclast precursors (pOC) and osteoclasts. Temporal changes in bone area (μ m 2 ) were also predicted under control and TGFβ inhibitor conditions. parameters and therefore, to compare computational model outputs to those obtained with the PAIII model, we used a scaling approach. Briefly, using tumor growth rates from the computational model (25 day intervals) and the animal model (2 day intervals) we calculated the slope of the line/derivative for each time point. Using this approach, we found that derivatives for each model are similar for pre and post treatment between day 8-10 for the animal model and day 75-100 for the in silico model. At these comparable points, the computational model accurately predicts the impact of TGFβ inhibition on cancer cell proliferation but differences in apoptosis were not evident in the computational model until later time points (Fig. 3d vs. 3g).
To further assess the predictive power of the computational model, we performed a TGFβ inhibitor post-treatment study. In silico findings suggest treatment of an established TGFβ ligand and receptor-expressing cancer (TRP) with a TGFβ inhibitor would not be of benefit. Our in vivo data confirm these predictions ( Supplementary Fig. S3).

TGFβ inhibition prevents prostate cancer induced osteolysis. Because of their role in the vicious
cycle, we next focused on osteoclast behavior over time in the control and TGFβ inhibitor treated groups. Dissection of the computational model outputs revealed the number of active bone resorbing osteoclasts were significantly lower in the TGFβ inhibitor versus control group at day 100 (Fig. 4a). Further, the rate of osteoclast maturation and apoptosis was significantly mitigated during this period, which is in agreement with known effects of TGFβ inhibition on osteoclast function (Fig. 4b,c) 28,33 . Analysis of the lytic component of our in vivo model identified that there was approximately 50% (p = 0.002) less tumor induced osteolysis in the TGFβ inhibitor group compared to the controls as determined by X-ray (Fig. 4d). Histochemical analysis using the osteoclast specific marker tartrate resistant acid phosphatase (TRAP) demonstrated that this reduction in osteolysis was due to significantly fewer osteoclasts in the TGFβ inhibitor treated group compared to the controls (Fig. 4e). These data suggest TGFβ inhibition reduces the extent of cancer-induced osteolysis over time in vivo by limiting osteoclast function and validates the pre-treatment HCA model outputs.
TGFβ inhibition has differential effects on normal and cancer induced bone formation. Based on published data demonstrating that TGFβ blockade increases bone formation 29 , we hypothesized that inhibition of the growth factor would significantly exacerbate prostate cancer induced osteogenesis. At day 100, the computational model does predict a small but significant increase in bone formation in the TGFβ inhibitor simulations despite a concomitant decrease in osteoblast proliferation (Fig. 5a-c and Supplementary Fig. S2). While microCT (μ CT) scans of tumor bearing bones did not reveal any differences between the groups (Fig. 5d), histomorphometry analysis does support the in silico predictions with increased bone formation in TGFβ inhibitor treated mice and a trend towards fewer bone rimming osteoblasts (Fig. 5e,f). In contrast, we observed a robust increase in bone formation in contralateral sham limbs derived from the TGFβ inhibitor group for all measurements ( Fig. 5d-f). Taken together, these in vivo data suggest that TGFβ inhibition does not greatly exacerbate prostate cancer induced osteogenesis compared to non pathological conditions and support, in part, HCA model outputs. It is also worth noting that over a longer period (Day 250), the computational model predicts TGFβ inhibition would ultimately result in decreased bone formation ( Fig. 5a and Supplementary Fig. S2), a finding that warrants further exploration in vivo.
Predicting TGFβ inhibitor efficacy on heterogeneous bone metastatic prostate cancer. While the majority of cancer cells in human specimens of bone metastatic prostate cancer produce both TGFβ ligand and receptors (TRP), we also noted the presence of cancer cells that either produced the ligand alone (TP), the receptor alone (TR) or neither (TN) (Fig. 2). This raises the question as to what the impact of TGFβ inhibition would be on these clonal populations. We have shown that, in our hands, the C4-2B cell line produces TGFβ but not the receptors (TP) (Fig. 2b-d). Therefore, we challenged the computational model to determine the impact of TGFβ inhibition on a homogenous TP bone metastatic cancer and found no effect on the growth of the cancer in pre-or-post treatment simulations (Equation 5, Supplementary Fig. S4). To test these results in vivo, C4-2B luciferase-expressing cells were inoculated into mice that were pre-or post-treated with TGFβ inhibitor or IgG control and we observed no difference in tumor growth between the groups (Supplementary Fig. S4).
Clinically, human samples of bone metastatic prostate cancer are heterogeneous for their usage of TGFβ (Fig. 2a). A major advantage of the HCA modeling approach is the ability to integrate multiple clonal phenotypes (g) The proliferative and apoptotic index in TGFβ i and control tissue sections were measured using pHistone H3 and cleaved caspase-3 (c-Caspase-3) (red) respectively as a ratio to total cell number (DAPI; blue). Scale bars, 50 μ m. Asterisks denote statistical significance (* p < 0.05; * * p < 0.005).
based on patient specific information. In this regard, analysis of a single patient specimen, identified heterogeneous expression and activity in the TGFβ signaling axis with the following clonal ratio noted: TP 1: TRP 231: TR 6: TN 4 (Fig. 6a,b). We seeded the HCA model with this ratio of clones and performed simulations for control, pre-and post-treatment conditions (n ≥ 24/group) to determine the impact of TGFβ inhibition on metastatic prostate cancer evolution (Equations 3-6). In silico data show the effects of pre-and post-TGFβ inhibitor treatment on clonal evolution for this patient over time and the responses of the surrounding tumor microenvironment ( Supplementary Fig. S5, Table S4, and Video 3). The model also predicts the dominance of the TRP clone over time under control conditions (Fig. 6c). However, application of the TGFβ inhibitor, especially in the pre-treatment groups at varying levels of efficacy, show how the cancer cell dynamics shift over time in favor of other clones, in particular for the TN population (Fig. 6d-f and Table S4). These data underscore the power of the HCA computational modeling approach in examining the temporal effects of targeted therapy on heterogeneous cancer cells and the surrounding microenvironment.

Discussion
Current treatment options for patients diagnosed with bone metastatic castrate resistant prostate cancer include second generation androgen deprivation therapies, radiation treatment, bisphosphonates/anti-RANKL therapies, alpharadin and/or chemotherapy 1 . While these treatments mitigate pain, pathological fracture and increase overall survival, the disease remains incurable with the median survival time being approximately 3-years subsequent to diagnosis. Increasing our knowledge of the mechanisms driving the disease can reveal novel therapeutic targets. To this end, dozens of molecular mechanisms that play important roles in tumor-bone interaction have been discovered 34 . However, dissecting how potential targeted therapies will work in the context of current treatment paradigms and their translation to the clinical setting presents a major challenge. Using TGFβ inhibition as an example, we demonstrated how a novel biologically driven computational HCA model can rapidly define temporal cancer-bone microenvironment responses to a given therapy. Further, the integrated approach provides insight into optimal therapeutic windows to apply a given inhibitor. For TGFβ inhibition, the HCA model indicates that application of the inhibitor in an adjuvant setting subsequent to the detection/treatment of aggressive prostate cancer would be most effective.
As with biological models, there are numerous mathematical approaches to study cancer progression such as branching and Moran processes based models, systems of ordinary and partial differential equations, and agent based models, each with varying strengths and weaknesses 6,15,35,36 . Data driven "top down" models "fit" existing clinical or experimental information to identify parameters that explain the behavior of the disease 37,38 . In the context of prostate cancer, many elegant models have used biological parameters such as prostate specific antigen (PSA) to predict time to progression, and to model the effects of intermittent androgen therapy 8 . The tightly regulated process of normal bone remodeling lends itself well to modeling how key factors such as RANKL and TGFβ control the behavior and activity of bone stromal cells over time [39][40][41][42][43] . By extension, perturbing this balanced ecosystem with an invasive species such as cancer can also be modeled. Our own group, as well as others, has been exploring how bone metastatic cells and skeletal malignancies such as myeloma interact with the bone microenvironment in order to progress 18,44,45 . In the current study, we used an agent-based approach that allows for the exploration of key cellular interactions over space and time in an unbiased hypothesis driven "bottom-up" manner. We have shown how experimentally-derived cellular parameters can be integrated into the rule sets and partial differential equations (PDEs) used to drive HCA models. This allows us to produce biologically 46 and clinically 47 testable hypotheses without making assumptions regarding population level dynamics. The result is a model capable of generating predictions that naturally emerge from the interactions between cells and their environment. Key to the robustness of the HCA model outputs, is the reliability of the biological parameters used to power the PDEs. In the current study, we used human parameters to generate the HCA model 21 . While we assume that individual bone stromal cell components have similar dimensions and lifespans between the species, human cancer cells typically grow much more slowly than human derived xenografts or animal models of the disease 21 . Therefore appropriate scaling between the models is an important consideration for the comparison of results. Further, not all of the HCA models predictions were correct. For example, we observed the in silico effects of TGFβ inhibition on osteoblast number was discordant with in vivo results (Fig. 5b vs. 5f), suggesting that the parameters/assumptions governing the effects of TGFβ inhibition on apoptosis require re-calibration based on the obtained biological data. This reiterative process allows for fine-tuning of the HCA.
The roles of TGFβ in skeletal development and malignancy have been well described 5,24 . However, enthusiasm for applying TGFβ inhibitors as a therapeutic strategy to treat metastatic bone disease including mCRPC is limited because of the pleiotropic and often opposing roles TGFβ plays in normal and cancer cell biology 5,24,25,48,49 . This complexity made TGFβ and the effects of TGFβ inhibition an ideal challenge with which to test the power of our HCA based computational model. Our results indicate the treatment of established active metastases with a TGFβ inhibitor, unless applied at > 99% efficacy, would have little or no impact on the progression of cancer cells regardless of their dependency on TGFβ . In contrast, the application of the inhibitor in a preventative or adjuvant manner would significantly control bone metastatic prostate cancer growth and osteogenesis. This result was predicated on the basis that the metastatic cancer cells have an active TGFβ signaling axis which we found to be the case in the majority of human bone metastatic prostate cancer specimens. The mathematical model was built on the assumption that TGFβ , RANKL and bone derived nutrients drive cell responses. Arguably then, interfering with TGFβ signaling in the model could potentially be self-fulfilling in predicting cancer-host behavior. However, despite the limited number of cytokines and growth factors included in the model, we validated many of the model predictions with independent in vivo experiments thus reinforcing the key roles for TGFβ and RANKL in the vicious cycle of tumor bone interaction. The majority of our in vivo results support the accuracy of the parameters and assumptions used to power the computational model. The in vivo results are consistent with other studies examining the role of TGFβ in bone metastatic cancers including prostate, breast and melanoma but, importantly, underscores the potential for computational modeling in predicting the efficacy of an applied targeted therapy [50][51][52] . Although our model is centered on the roles of TGFβ and RANKL, it is important to note that the circuitry of the HCA model can easily be expanded to include other cell types (immune cells or cancer-associated fibroblasts) and molecules (e.g. PTHrP, BMPs) and their roles/effects subsequently explored. While the computational model was primarily used to examine TGFβ inhibition in bone metastatic prostate cancer, the outputs generated also revealed new insights into TGFβ biology. Firstly, the model highlights the critical role for TGFβ in the cyclical dynamics of populations such as osteoclasts and confirms the importance of the osteoclast in initializing the vicious cycle. For example, we noted osteoclast infiltration often precedes a period of cancer growth and osteogenesis (Video S2), a finding that could be informative for the timing of anti-resorptive therapies as we have previously shown 18 . TGFβ inhibition alters these dynamics and at specific time-points can reverse the trends between the treated and control groups so that osteoclast numbers in the control simulations may in fact be lower than those in the TGFβ inhibitor treated group (Fig. 3A). This suggests that arbitrarily selected time-points in pre-clinical in vivo animal studies may not accurately reflect how applied therapeutics are impacting cancer-bone interaction over time. Secondly, although TGFβ inhibition promotes robust bone formation in normal non-pathological situations, the model and our in vivo results confirm counter-intuitively that, TGFβ inhibition does not greatly exacerbate prostate cancer induced osteogenesis.
Currently, mCRPC inter-patient heterogeneity is a major clinical challenge. Integrating biological and computational modeling offers a unique opportunity to study how cancer evolves and reacts to changing microenvironments and applied targeted therapies. Patient derived xenograft (PDX) samples are being used to design precision treatment strategies and integrating the biological parameters derived from these specimens into computational models could prove to be a synergistic way to tackle the complexity of heterogeneity in individual patients 19 . Again, using patient derived TGFβ signaling axis information as an example, we demonstrated how clonal variation and evolution in response to applied TGFβ inhibitors could be incorporated into the HCA model (Fig. 6). In these studies, we assumed the TGFβ dependent growth rates of clonal cancer cells. For personalization of the HCA and its clinical application, individual patient specimens would have to be isolated and growth rates of various clones examined in ex vivo assays. Current advances in single cell ex vivo analyses support the feasibility of such an approach 53,54 . The HCA model also allows for the optimization of inhibitor dosing and timing that in turn could be used to generate an adaptive therapy strategy to prevent the outgrowth of resistant sub-populations 20,55 .
In conclusion, we have developed a novel and unbiased computational HCA model that allows for the dynamic multi-scale understanding of how metastatic prostate cancer cells evolve and interact with the surrounding bone microenvironment. We used this model to predict the efficacy and response of bone metastatic prostate cancer to targeted therapies such as, TGFβ inhibitors. Further, our integrated computational and biological approach allowed for the dissection of how TGFβ inhibition simultaneously affects osteoblast, osteoclast and cancer cell behavior over time. The HCA model constitutes a platform for discovery that can readily be expanded to incorporate additional cellular and molecular circuitry. This will ultimately yield a clinical tool that will aid the medical oncologist in designing curative strategies for heterogeneous bone metastatic prostate cancer. Most importantly, this modeling approach can be applied to the development of new therapeutic strategies across a broad spectrum of human malignancies.

Methods
An HCA Model of bone metastatic prostate cancer for therapy optimization. In the computational model, each cell type responds to TGFβ levels in an either directly proportional (1 + Log(TGFβ)) or inversely proportional (− 1Log(TGFβ)) manner. TGFβ inhibition in the HCA is achieved by controlling TGFβ bioavailability (0-99% inhibition at a constant level). Computational models were seeded with homogeneous or heterogeneous prostate cancer cells that expressed the TGFβ receptor and ligand (TRP), the ligand alone (TP), the receptor alone (TR) or were negative for both (TN). For the HCA model, we consider cellular intrinsic behaviors and the impact of TGFβ on these behaviors. We include 8 different cell types: resident and active cells of the bone stroma (mesenchymal stromal cells (MSCs), precursor osteoblasts (pOB), adult osteoblasts (aOB), precursor osteoclasts (pOC), adult osteoclasts (aOC)) and prostate cancer cells of varying TGFβ responsiveness (TRP, TGFβ ligand and receptor-producing; TR, receptor-producing; TP, ligand-producing and; TN, negative for receptor and ligand expression). We considered interactions between all cell types and the impact of these interactions on the tumor-bone microenvironment. Empirical, experimental and theoretical parameters were used to fuel equations (Table S1).
Precursor Osteoblasts and Osteoclasts. Based on empirical data and literature, we assume a maximum rate of pOB division as 36 hours based on ATCC specifications for MC3T3-E1 cells. The rate of pOB division is inversely affected by TGFβ 28 and we assume that the effect is logarithmic. If TGFβ is at saturation levels (> 10 ng/ml), then the pOB division rate tends to zero. This assumption is based on our findings of the effect of TGFβ on osteoblast precursors (Supplemental Fig. S1). If there is no TGFβ present, then the maximum rate is considered. The rate of division has subsequent effects on the number of mature bone generating osteoblasts. By the same token, the fusion rate of pOC is also affected inversely proportional to the availability of TGFβ . These behaviors can be described by: precursor where 0 < TGFβ < 1 and Div precursor is substituted in the case of precursor osteoblasts, Div pOB is the maximum rate of pOB division, and in the case of precursor osteoclasts, Div pOC , the maximum rate of pOC fusion in the absence of TGFβ . This ensures that, when bone is being resorbed and TGFβ is being made bioavailable, osteoclastogenesis is limited. Once fused, the probability of aOC survival also depends proportionally on TGFβ . This is calculated by: aOC where 0 < TGFβ < 1 and Surv aOC is the maximum percentage of survival for aOCs (100% when TGFβ is as its saturation level). If the levels of TGFβ are below the saturation level, the probability of death for the aOC increases 56,57 .
Bone metastatic prostate cancer cells. Prostate cancer cells were explicitly defined as being dependent on TGFβ and bone derived nutrients (BDN) for their division. The probability of division was estimated as being proportional to the inverse logarithm of the available BDN. If there are no nutrients, there is zero division. If there is maximum nutrient saturation, the division rate is at its maximum. Based on empirical data obtained with TRP cells such as PAIII, we assume that TRP cells have a maximum division rate of once every 1.5 days and a lifespan without contact with BDN of 14 days based on low serum (2%) soft agar assays ( Supplementary Fig. S4). TRP division in response to bone derived nutrients was modeled as follows: TRP where 0 < TGFβ < 1 and Div TRP is the maximum rate of TRP division, in the saturation level of bone derived nutrients. Under TGFβ inhibition, the maximum division rate is assumed to be reduced for TRP to once every 2 days based on empirical observations in our laboratory (data not shown). We assume that TR has a maximum division rate of 1.75 days, but a lifespan without nutrients of 10 days, calculated as not having a cost for producing TGFβ ligand but benefiting from the presence TGFβ in the bone microenvironment. TR cell division depends directly on bone-derived nutrients in the same manner: TR where 0 < BDN < 1 and Div TR is the maximum rate of TR division, in the saturation level of bone derived nutrients. Under TGFβ inhibition, the maximum division rate is reduced for TR to once every two days. We assume that TP has a maximum division rate of once 1.75 days, and a lifespan without nutrients of 10 days, calculated as having a cost for producing TGFβ ligand but not benefiting from TGFβ contained within the bone derived nutrients. TP cell division depends directly on bone-derived nutrients in the same manner: Finally, we assume TN has a maximum division rate of once every 2 days and a lifespan without nutrients of 12 days, calculated from having no cost of producing TGFβ and no benefit from TGFβ signaling. TN cell division depends directly on bone-derived nutrients through: Cell Culture and Patient Specimens. Luciferase-expressing PAIII, C42B, and PC3 prostate cancer cell lines were cultured in complete Dulbecco's Modified medium supplemented with 10% fetal bovine serum [58][59][60] .
All cell lines were periodically tested for mycoplasma (#CUL001B, R&D Systems) and short tandem repeat (STR) verified at the Moffitt Clinical Translational Research Core. De-identified tissue sections of bone metastatic prostate cancer were obtained from the Moffitt tissue archives (MCC 50086).
In vivo experiments. All animal experiments were performed with IACUC approval (R1283) and were conducted in accordance with the guidelines set forth in the Guidelines for the Care and Use of Laboratory Animals published by the National Institutes of Health. Pre-Treatment Studies: 6-week old male SCID Beige mice were injected intraperitoneally with either TGFβ inhibitor (1D11; 10 μ g/ml; n = 10/group) or isotype control, 13C4 (10 μ g/ml; n = 8/group), a kind gift from Scott Lonning and Patrick Finn at Genzyme. Subsequently, luciferase expressing PAIII or C4-2B cell lines (5 × 10 4 or 1 × 10 5 respectively in 20 μ l of saline) were intratibially injected either one-day or one week after TGFβ pre-treatment 61,62 . Contra-lateral limbs were injected with saline and served as a positive control. Mice received TGFβ inhibitor or IgG control injections every three days (PAIII model) or weekly (C4-2B). Bioluminescence was measured longitudinally as a correlate of tumor growth (IVIS ™ Perkin Elmer). For Post-treatment studies: mice were inoculated as described (n = 7/group), randomized and treated upon the detection of bioluminescent signal. Mice that showed tumor growth outside of the bone compartment were excluded from all analyses.
A detailed description of cell assays, histological and bone morphology analyses can be found in supplementary methods. Statistical analyses were performed with GraphPad Prism and all graphs display error bars that are SEM.