\documentclass[11pt]{article} \setlength{\topmargin}{-0.3in} \renewcommand{\baselinestretch}{1.5} \renewcommand{\arraystretch}{0.6} \setlength{\topskip}{0.1in} \setlength{\textheight}{9.2in} \setlength{\oddsidemargin}{0.1in} \setlength{\evensidemargin}{0.1in} \setlength{\textwidth}{6.75in} \usepackage{amsmath} \usepackage{float,amsthm} \usepackage{amsfonts} \usepackage{siunitx} \usepackage[colorlinks=true,citecolor=blue]{hyperref} \usepackage[mathlines,displaymath]{lineno} %\runninglinenumbers \usepackage[title,toc]{appendix} % 'title' adds "Appendix" before title, 'toc' adds to TOC \usepackage[english]{babel} \usepackage{float} \usepackage{epsfig} \usepackage{epstopdf} \usepackage[toc]{appendix} \usepackage{graphicx} \newcommand{\R}{\mathbb R} \newcommand{\C}{\mathbb C} \newtheorem{corollary}{Corollary} \newcommand{\ds}{\displaystyle} \newtheorem{theorem}{Theorem}[section] \newtheorem{prop}{Proposition}[section] \newtheorem{definition}{Definition}[section] \newtheorem{remark}{Remark}[section] \newtheorem{lemma}{Lemma}[section] %\usepackage{mathtools} \numberwithin{equation}{section} \begin{document} \title{\bf Computational Models of Lung Cancer Cell Dynamics and Immune System Interplay} \vspace{6in} %\author{{\normalsize Salamida Daudi$^{1}$\footnote{Correspondence email.: daudisalamida81@gmail.com; mhelikumi@yahoo.co.uk; steadymushaya@gmail.com},~~~Mlyashimbi Helikumi$^{2},$ ~Steady Mushayabasa$\,^{3}$} \\ % {\scriptsize \it $^{1}\,$ Department of Mathematics, Humanities and Social Science, National Institute of Transport, Dar-es-Salaam, Tanzania}\\ % {\scriptsize \it $^{2}\,$ Mbeya University of Science and Technology, Department of Mathematics and Statistics,}\\ %{\scriptsize \it College of Science and Technical Education, P.O. Box 131, Mbeya, Tanzania,}\\ %{\scriptsize \it $^{3}\,$University of Zimbabwe, Department of Mathematics \& Computational Sciences, Harare, Zimbabwe} %\\ %{\scriptsize P.O. Box MP 167 Mount Pleasant, Harare, Zimbabwe} %\vspace{0.1in}} \date{} \maketitle \begin{center} \noindent{\bf Willie K Chidaushe . Thomas Musora. Steady Mushayabasa} \end{center} \hrule \begin{center} \begin{abstract} \noindent Although lung cancer is one of the types of cancer leading to death worldwide, many investigations have been performed. Mathematical models were formulated to explain the mechanisms that are crucial to predict the behavior of cancer cell growth in the lungs, invasion of cancer to tissue, and the effects of anticancer agents after administration together with response of the immune system. The literature on mathematical modeling of lung cancer dynamics involving invasion and migration is abundant. Mathematical models to simulate the growth rate of the cancer cells derived from both deterministic and stochastic considerations. The anti cancer model formulated and analyzed from the exponential basic equation explaining the growth of tumor cells was considered the basic beginning of the model. Since the sensitivity analyses normally suggest that the anti cancer agent reduces tumor mortality rate and drug decay rate will contribute to the determination of treatment outcomes.\\ \vspace{0.7cm} \noindent{\bf Keywords}: Lung cancer, Mathematical model, Immune system interplay , Treatment, Drug resistance\\ %\noindent{\bf MSC codes.:} 35B40, 35K57, 35Q92, 92D30 \end{abstract} \end{center} \hrule \section{Introduction} Cancer is the deadliest and complicated disease of our time. This illness is caused by the uncontrolled growth of abnormal or mutated cells in the body, and cancer cells are known for their capability to grow rapidly, divide and proliferate uncontrollably \cite{ozkose2021fractional}. The cancer development is due to tumor cell growth which is caused when normal cells become cancerous when a series of mutations leads cells to continue to grow and divide out of control.Solid tumors contain multiple mutations also known as variants, which represent a change from the original. Cancer cells are different from normal cells and they have the ability to invade nearby tissues and spread to distant regions of the body. They appear through a series of genetic and environment-induced changes (epigenetic), which may be inherited or, more often, caused by carcinogens or cancer-causing substances in the environment. The proposed research will employ tumor growth in thorax or lungs area of human body as main case study. Lung cancer is a kind of cancer that starts as a growth of tumor cells in the lungs \cite{tao2019epidemiology}. The lungs are two spongy organs in the chest that control breathing and are located in the thoracic cavity resting on a muscle called the diaphragm. People who smoke have the greatest risk of lung cancer. The risk of lung cancer increases with the length of time and number of cigarettes smoked. Quitting smoking, even after smoking for many years, significantly lowers the chances of developing lung cancer. Lung cancer can also affect people who have never smoked. \\ Normal cells in the thorax grow during development stages like in other body parts, such as during childhood, or to repair injured tissue, while cancer cells continue to grow even when more cells are not needed. Normal lung cells are responsible for the crucial function of gas exchange, facilitating the transfer of oxygen into the bloodstream and carbon dioxide out of the body. They also play a role in maintaining the acid-base balance of the body and are involved in the body's immune response. Cancer cells grow continually without a limit and also fail to listen to signals that tell them to stop growing, which even leads them to commit suicide. These cells can be affected by various factors, including exposure to toxins, which can lead to damage and potentially cancerous changes \cite{hsia2016lung}. This cells differ that the others they have the ability to invade nearby cells but normal cells respond to signals from other cells which tell them that they have reached a boundary of population growth. Tumor cells extend into nearby tissues, often with finger like projections. Cancer has many causes both inside and outside the body contribute to the development of cancer \cite{boffetta2003contribution}. Environmental factors are factors outside the body which includes cigarette smoking, excessive alcohol consumption, poor diet, lack of exercise, excessive sunlight exposure, and sexual behavior that increases exposure to certain viruses.\\ The development and healthy life of a human being requires the cooperation of more than ten million cells for the good of organism \cite{venter2004century}. This co-operation is maintained by signals and cellular checkpoints which determine whether cell divide, die or differentiate. Cancer represents the collapse of this cooperation this results in a uncontrolled growth of cells within the body which eventually leads to the death of the organism \cite{tripathy2020cancer}. The uncontrolled growth of cells is the results to alterations or mutations in the genetic materials. This abnormal cell behavior disrupts normal tissue function and can lead to various complications and, ultimately, death. More precisely, the emergence of cancer may require the accumulation of multiple mutations which allow cell to break out of the regulatory networks which ensure cooperation which refers to as a multi-stage carcinogenesis. Cancer is the second-leading cause of death in the world. But survival rates are improving for many types of cancer, thanks to improvements in cancer screening, treatment and prevention.\\ Through many researches made for many years to prevent lung cancer growth and recurrence, many anti cancer agent were discovered which are also known as anti cancer drugs or cancer therapeutics substances that inhabits or stops the growth and proliferation of tumor cells \cite{devita2012two}. They are several types of anti cancer agents which can be used for treatment to reduce the tumor development or at some point will end tumor development at earliest stage of cancer development which are chemotherapy agent, targeted therapies, hormonal therapies, immunotherapy, radiotherapy enhancing agents, gene therapy agents and plant derived anti cancer agent. These anti cancer agents can be combined for treatment to produce various effects, on the positive side of effects some can enhance anti tumor efficacy improving overall survival while others can causes negative effects like increase of toxicity and can cause risk of adverse interactions.\\ Lung cancer is a type of cancer that starts also when abnormal cells grow in an uncontrolled way in the lungs \cite{pandi2016brief}. It is a serious health issue that can cause serious harm and death. Cancer diseases are leading to death worldwide, with severe socioeconomic repercussions \cite{vineis2014global}. To better understand these repercussions, they is need to investigate the tumor growth kinetics and the effects of co-administered anti cancer agents describing the limited growth in number of tumor cells, using mathematical models. A predictive model of tumor growth kinetics can lead to early detection of cancer development as well as reducing financial burden for both the patient and society. \subsection{Background} The first task faced by medicine was discovered clinically and macroscopically in different forms of cancer and to unite them in some common notion; more recently some knowledge as to its microscopic morphology, its physiology and causative factors has been acquired and some new therapeutic agents have been discovered. Lung cancer was first described by doctors in the mid-19th century. In the early 20th century it was considered relatively rare, but by the end of the century it was the leading cause of cancer-related death among men in more than 25 developed countries \cite{rosenberg2007erwin}. Lung cancer,is also called lung carcinoma, is a malignant tumor that originates in the tissues of the lungs. In the 21st century, lung cancer has emerged as the leading cause of cancer deaths worldwide. By 2012, it had outpaced breast cancer as the leading cause of cancer death among women in developed countries. The rapid increase in the worldwide prevalence of lung cancer was mainly attributed to the increased use of cigarettes after World War I, although increases in environmental air pollution were suspected to have been a contributing factor as well.\\ Lung cancer was rare before the advent of cigarette smoking. Surgeon Alton Ochsner recalled that as a Washington University medical student in 1919, his entire medical school class was summoned to witness an autopsy of a man who had died from lung cancer, and told they may never see such a case again\cite{spiro2005one}. In Isaac Adler's 1912 Primary Malignant Growths of the Lungs and Bronchi, he called lung cancer "among the rarest forms of disease"; Adler tabulated the 374 cases of lung cancer that had been published to that time, concluding the disease was increasing in incidence.By the $1920s$, several theories had been put forward linking the increase in lung cancer to various chemical exposures that had increased including tobacco smoke, asphalt dust, industrial air pollution, and poisonous gasses from World War I. Lung cancer today is primarily caused by the inhalation of smoke from cigarettes, which is also why the disease was quite rare prior to the 20th century. Lung cancer was not even recognized medically until the 18th century, and as recently as 1900 only about 140 cases were known in the published medical literature. The findings of primary lung tumors in autopsied bodies of German research clinics rose dramatically in the second half of the nineteenth century and even more dramatically in the first decade of the twentieth century \cite{timmermann2014lung}. Isaac Adler summarized this evidence in 1912, in the world's first monograph on lung cancer, noting that the incidence of malignant neoplasms of the lung seemed to show ‘a decided increase’ \cite{proctor2012history}. Adler mentioned the ‘abuse of tobacco and alcohol’ as one possible cause, while also commenting that the subject was ‘not yet ready for final judgment’.\\ By late $1953$, the tobacco industry faced a crisis of cataclysmic proportions. Smoking had been categorically linked to the dramatic rise of lung cancer. Although health concerns about smoking had been raised for decades, by the early $1950s$ there was a powerful expansion and consolidation of scientific methods and findings that demonstrated that smoking caused lung disease as well as other serious respiratory and cardiac diseases, leading to death \cite{proctor2012history}. These findings appeared in major, peer-reviewed medical journals as well as throughout the general media. At the same time, internal research at the major tobacco companies supported the link between tobacco and lung cancer; though these results were kept secret from the public \cite{proctor2012history, bates2004tobacco}. The connection of lung cancer with radon gas was first recognized among miners in Germany's Ore Mountains. As early as 1500, miners were noted to develop a deadly disease called "mountain sickness" ("Bergkrankheit"), identified as lung cancer by the late 19th century \cite{mc2012historical}. In the $1950s$ radon and its breakdown products became established as causes of lung cancer in miners. Based largely on studies of miners, the International Agency for Research on Cancer classified radon as "carcinogenic to humans" in 1988. In 1956, a study revealed radon in Swedish residences. During the following decades, high concentrations of radon were found in residences throughout the world; By $1980s$, many countries had established national radon programs to catalog and mitigate residential radon. \subsection{Existing Treatment Options of Lung Cancer} They are the primary treatment options for cancer which includes surgery, radiation therapy, and systemic therapy (including chemotherapy, targeted therapy, immunotherapy, and hormone therapy). Recently, many pathways involved in cancer therapy progression and how they can be targeted has improved dramatically, with combinatorial strategies, involving multiple targeted therapies or “traditional” chemotherapeutic, such as the taxanes and platinum compounds, being found to have a synergistic effect \cite{debela2021new}. New approaches, such as drugs, biological molecules, and immune-mediated therapies, are being used for treatment even if the excepted therapy level has not reached that resists the mortality rate and decreases the prolonged survival time for metastatic cancer. These treatments can be used alone or in combination, depending on the type and stage of the cancer, as well as the patient’s overall health. The types of treatment that you receive will depend on the type of cancer you have and how advanced it is but for lung cancer the treatment depends with type. Lung cancer is categorized into two type which are non small cell lung cancer and small lung cancer each with different treatment options.\\ Non-small cell lung cancer (NSCLC) occurs when normal cells in your lungs change and grow out of control. NSCLC grows slowly compared to small cell lung cancer. But it can spread to other parts of your body before you develop noticeable symptoms. Early detection and treatment are key. Early stage of non small lung cancer may be treated with surgery followed by chemotherapy (adjuvant chemotherapy). If surgery is not possible, appropriate or acceptable, a combination of radiation treatment and chemotherapy may be recommended \cite{burdett1996adjuvant}. Some types of advanced non-small cell lung cancer may respond to immunotherapy treatment. Targeted therapy can be effective for people with advanced cancer who have a specific gene change. However small cell lung cancer is a fast growing type of lung cancer, characterized by small, oval-shaped cells when viewed under a microscope. It is an aggressive form of cancer, meaning it tends to spread rapidly to other parts of the body \cite{cooper2006small}. Immunotherapies have emerged as a highly effective treatment option for hematological malignancies, with Chimeric Antigen Receptor (CAR)T-cell therapy being the most successful in use today. In these innovative therapies, first approved by the FDA (Food and Drug Administration) only in 2017, T-cells are extracted from the patient’s blood and undergo genetic engineering within a laboratory to introduce a chimeric antigen receptor (CAR) tailored to cancer cells \cite{sabir2025mathematical}. Subsequently, these modified T-cells are cultured and expanded, generating a robust population. Once infused back into the patient, the CAR T-cells recognize and bind to cancer cells by targeting specific proteins on their surface. This binding activates the CAR T-cells, initiating a robust immune response characterized by the release of cytotoxic substances, ultimately leading to the destruction of cancer cells. \section{ Models and Methods} The modeling of infectious diseases is a tool that has been used to study the mechanisms by which diseases spread, to predict the future course of an outbreak and to evaluate strategies to control an epidemic. mathematical modeling of spread of disease was carried out in 1760 by Daniel Bernoulli. Trained as a physician, Bernoulli created a mathematical model to defend the practice of inoculating against smallpox \cite{hethcote2000mathematics}. The calculations from this model showed that universal inoculation against smallpox would increase the life expectancy from 26 years 7 months to 29 years 9 months. Daniel Bernoulli's work preceded the modern understanding of germ theory. Mathematical models play a crucial role in lung cancer research by providing a framework to understand, predict, and potentially control cancer growth and treatment response. The formulated models can simulate tumor behavior, optimize treatment strategies, and offer insights into the complex dynamics of cancer progression.\\ Mathematical modeling contributes to cancer research by helping to elucidate mechanisms and by providing quantitative predictions that can be validated \cite{altrock2015mathematics}. Various mathematical model structures have been used to characterize the tumor dynamics and drug resistance evolution for solid tumors in lungs. Tumor proliferation, regression due to treatment, heterogeneity, and treatment resistance are key elements that are commonly considered in those models. Many mathematical and computational models have been developed to simulate tumor growth lungs and drug response. Therefore, the development of mathematical models capable of quantitatively evaluating synergism in combination drug therapy is desirable. \subsection{Mathematical Model} We constructed an ordinary differential equations (ODEs) model to understand the implications of combining drugs in limiting lung cancer resistance. The proposed model considers two tumors, the tumor cells $C(t)$ described as abnormal cells or premalignant cells that can develop into a clones of cancer cells and then the tumor cell clones $C_{N}(t)$ which are genetically identical groups of tumor cells that have all descended from a single original cancer cell. These clones arise through mutations and compete with each other within a tumor for resources \cite{fialkow1976clonal}. Additionally, When a lung tumor develops, the expression of estrogen receptors (ERs) denoted by $E(t)$ is altered compared to healthy lung tissue, most notably by the prominent expression of estrogen receptor beta (ER$\beta$), which is associated with tumor growth and metastasis. A tumor can cause an obstruction in the airway also bronchial obstruction which causes lung tissues to breaks down, leading to the accumulation of endogenous lipids (fats produced by the body) and lipid-laden macrophages, resulting in endogenous lipoid pneumonia denoted by $F$ \cite{alvarado2016metabolic}. Thus, to model the immune response let the variable $T_{c}(t)$ denote the population of effector T cells like Innate Lymphoid Cells (CTL) left by the naive T cells which then proliferate and differentiate into armed effector T cells to find and eliminate the source of the antigen \cite{xiong2025mathematical} . Antigen generation starts when cancer cells in the lungs accumulate genetic mutations, leading to the production of abnormal or over-expressed proteins not found in normal healthy cells. These antigens activate the free endogenous antigenic peptides denoted by $P$. Macrophages in lung cancer known as Tumor-Associated Macrophages (TAMs) with two phenotypes $M_{1}$ and $M_{2}$ denoted with $u_{M_{1}}$ and $u_{M_{2}}$ respectively are to be included in the model because they are the key player in promoting tumor growth, invasion, metastasis, and treatment resistance , primarily by adopting an $M_{2}$-like phenotype that supports angiogenesis, immune evasion, and matrix remodeling, though they can also have anti-tumor functions \cite{eftimie2021mathematical}. Additional assumptions are as follows: \begin{description} \item (i) The tumor cell population follows the logistic growth law \cite{xiong2025mathematical,wang2025optimal}. \item (ii) The macrophages phenotypes with highly plastic immune cells follows a logistic growth \cite{eftimie2021mathematical}. \item (iii) Effector T cells clears the damaged DNA \cite{Akman,Doisneau2023}. \item (iv) Production of estrogen which binds receptors promoting the growth and survival of cancer cells in lungs \cite{hsu2017estrogen,siegfried2009estrogen}. \item (v) The growth of lipoid pneumonia volume satisfies the logistic growth model \cite{Akman,Carrillo}. \item (vi) All tumor cells (tumor cells from damaged DNA and tumor cell clones) consume fat as energy resource \cite{Akman,Carrillo}. \item (vii) Activation of effector T cells by the dendritic cells and clearance of the effector T cells by the tumor cell clones \cite{xiong2025mathematical,ni1997role}. \item (vii) Immature dendritic cells becomes mature by capturing free antigens in tissues and antigens on the surface of living tumor cells. \cite{xiong2025mathematical,wang2025optimal}. \item (ix) The free antigenic peptides are derived from lysed tumor cells\cite{xiong2025mathematical}. \end{description} Based on these assumptions, the mathematical representation of the dynamics of cancer cells and immune response is described by a system of nonlinear ordinary differential equations (\ref{eq13})-(\ref{eq21}): \newpage \begin{eqnarray} \ds \frac{dC}{dt}&=&\ds \underbrace{\gamma _{C}C \left( {\begin{array}{c} 1-\frac{c}{K_{C}}\end{array} } \right)(1+r_{m_{2}}u_{M_{2}})}_\textrm{Logistic growth term for tumor cells}-\underbrace{\beta _{C}\frac{T^{2}_{c}C}{T^{2}_{c}+M_{T_{c}}}}_\textrm{$\left. \begin{tabular}{ll} Tumor cells killed by \\ effector T cells \end{tabular} \right. $ }\cr && - \underbrace{\beta _D\frac{D_{T}C}{D_{T}+M_{D}}}_\textrm{$ \left. \begin{tabular}{ll} Clearance of tumor cells by\\ mature dendritic cell \end{tabular} \right. $ }+\underbrace{\frac{\eta _{C}N_{0}C}{n}}_\textrm{$ \left. \begin{tabular}{ll} Transmission of damaged\\ DNA to tumor infected cells \end{tabular} \right. $ }-\underbrace{d_{C}u_{M_{1}}C}_\textrm{Natural decay}~~~~~~~~~~~~~~~~\label{eq13}\\[5pt] \ds \frac{dC_{N}}{dt}&=&\ds \underbrace{\gamma _{C_{N}}C_{N}\left( {\begin{array}{c} 1-\frac{C_{N}}{K_{C_{N}}}\end{array} } \right)(1+r_{M_{2}}u_{M_{2}})}_\textrm{Logistic growth term for tumor cell clones}-\underbrace{\beta _{C_{N}}\frac{C_{N}T^{2}_{c}}{T^{2}_{c}+M_{T_{c}}}}_\textrm{$\left. \begin{tabular}{ll} death of tumor cells by \\ effector T cells \end{tabular} \right. $ }\cr &&-\underbrace{\beta _D\frac{C_{N}D_{T}}{D_{T}+M_{D}}}_\textrm{$ \left. \begin{tabular}{ll} clearance of tumor\\ cells by dendritic cells \end{tabular} \right. $ }+\underbrace{\frac{\eta _{C_{N}}N_{0}C_{N}}{n}}_\textrm{$ \left. \begin{tabular}{ll} Transmission of damaged DNA\\ to tumor infected cell clones \end{tabular} \right. $ }-\underbrace{(e_{C_{N}}+d_{C_{N}})u_{M_{1}}C_{N}}_\textrm{Natural decay and death}\label{eq14}\\[5pt] \ds \frac{du_{M_{1}}}{dt}&=&\underbrace{\gamma _{u_{M_{1}}}u_{M_{1}}\left( {\begin{array}{c} 1-\frac{u_{M_{1}}+u_{M_{2}}}{K_m}\end{array} } \right)}_\textrm{$ \left. \begin{tabular}{ll} Logistic growth term for \\ macrophages with phenotype $M_{1}$ \end{tabular} \right. $ }-\underbrace{\alpha _{m_{1}}u_{M_{1}}\frac{C+C_{N}}{C+C_{N}+M_{u_{M_{1}}}}}_\textrm{ $ \left. \begin{tabular}{ll} Polarization of macrophages\\ with phenotype $M_{1}$ by tumor \end{tabular} \right. $ }-\underbrace{d_{m_{1}}u_{M_{1}}}_\textrm{Natural death}\label{eq15}\\[5pt] \ds \frac{du_{M_{2}}}{dt}&=&\underbrace{\gamma _{u_{M_{2}}}u_{M_{2}}\left( {\begin{array}{c} 1-\frac{u_{M_{1}}+u_{M_{2}}}{K_m}\end{array} } \right)}_\textrm{$ \left. \begin{tabular}{ll} Logistic growth term for \\ macrophages with phenotype $M_{2}$ \end{tabular} \right. $}-\underbrace{d_{m_{2}}u_{M_{2}}}_\textrm{Natural death}-\underbrace{\alpha _{m_{2}}u_{M_{2}}}_\textrm{ $ \left. \begin{tabular}{ll} Re-polarization of the $M_{2}$\\ macrophages \end{tabular} \right. $}\label{eq16}\\[5pt] \ds \frac{dE}{dt}&=&\underbrace{\lambda _{E} E_{0}\left( {\begin{array}{c}\frac{u_{M_{1}}+u_{M_{2}}}{u_{M_{1}}+u_{M_{2}}+M_{E}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{E}}{E+M_{E}}\end{array} } \right)}_\textrm{Activation of estrogen by macrophages}+\underbrace{(1-u)r_{E}F}_\textrm{$ \left. \begin{tabular}{ll} Estrogen receptors \\ production \end{tabular} \right. $}-\underbrace{\mu E}_\textrm{Degradation of ERs}\label{eq17}\\[5pt] \ds \frac{dF}{dt}&=&\underbrace{r_{F}F\left( {\begin{array}{c}1-\frac{F}{F_{Max}}\end{array} } \right)}_\textrm{Logistic growth of fats volume}-\underbrace{\alpha (C+C_{N})}_\textrm{Energy consumption}\label{eq18}\\[5pt] \ds \frac{dP}{dt}&=&\underbrace{k_{1}\beta _{P}\left( {\begin{array}{c}\frac{T^{2}_{c}}{T^{2}_{c}+M_{T_{c}}}\end{array} } \right)(C+C_{N})}_\textrm{ necrotic tumors induced}-\underbrace{d_{P} P}_\textrm{Natural decay}\label{eq19}\\[5pt] \ds \frac{dD_{T}}{dt}&=&\underbrace{\lambda_{D}D_{T0}\left( {\begin{array}{c}\frac{P+k_{2}(C+C_{N})}{P+k_{2}(C+C_{N})+M_{P}}\end{array} } \right)}_\textrm{ Immature dendritic cells becomes mature}-\underbrace{d_{D} D_{T}}_\textrm{death of mature dendritic cells}\label{eq20}\\[5pt] \ds \frac{dT_{c}}{dt}&= &\underbrace{\lambda_{T_{c}}T_{0}\left( {\begin{array}{c}\frac{D_{T}}{D_{T}+M_{D}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{T_{c}}}{T_{c}+M_{T_{c}}}\end{array} } \right)}_\textrm{ Activation of effector T cells}-\underbrace{\frac{\beta _{T_{c}}C_{N}T_{c}}{C_{N}+M_{C_{N}}}}_\textrm{Clearance by tumor cell clones}-\underbrace{d_{T_{c}} T_{c}}_\textrm{death of effector cells}\label{eq21} \end{eqnarray} with initial conditions: \begin{equation} C(0) =C_{0},\quad C_{N}(0) = C_{N0},\quad u_{M_{1}}(0) =u_{M_{1}0},\quad u_{M_{2}}(0) =u_{M_{2}0},\quad E(0) = E_{0},\quad \end{equation} And \begin{equation} F(0) = F_{0},\quad P(0) =P_{0},\quad D(0) =D_{0},\quad T_{c}(0) = T_{c0}. \label{eq2} \end{equation} In Eq (\ref{eq13}) were the tumor cells population follows a logistic growth with $\gamma _{C}$ ($0 \leq \gamma _{C}\leq 1$) as the intrinsic rate, $K_{C}$ as the carrying capacity, $r_{m_{2}}$ the rate at which the macrophages with phenotype $M_{2}$ contributes to proliferation of tumor cells and $r_{m_{12}}$ rate at which the mixed phenotypes ($M_{1}$ and $M_{2}$) with an assumed constant value denoted with $u_{M_{12}0}\approx 0$ meaning the ability of the immune system to transition between inflammation and healing is crippled by tumor in lungs. The tumor cells killed by effector T cells are assumed to satisfy the Hill equation with Hill coefficient $n=2$, that is, $\beta_{C}\frac{T^{2}_{c}C}{T^{2}_{c}+M_{T_{c}}}$ were $\beta _C$ representing the maximum killing rate followed by tumor cleaning by dendritic cells also modeled from the Michaelis-Menten term with the parameter of tumor cell cleaning $\beta _{D}$ \cite{wang2025optimal,xiong2025mathematical}. The normal cells are in the susceptible class denoted with a constant $N_{0}$ when they accumulate genetic mutations (errors in DNA) that disable growth controls (like tumor suppressor genes), lose normal cell death signals, and fail to stop dividing and often aided by a supportive tumor microenvironment with probabilities denoted by $\eta _{C}$ and $\eta _{C_{N}}$, allowing them to multiply uncontrollably and form masses. \\ Tumor cells turn into tumor clones through a process called clonal evolution when a single cell acquires a mutation which allows it to divide and create a clone of itself. Eq (\ref{eq14}) represents the tumor cell clones $C_{N}$ with a logistic growth were $\gamma _{C_{N}}$ ($0 \leq \gamma _{C_{N}}\leq 1$) is the intrinsic rate and $K_{C_{N}}$ is the carrying capacity. The tumor cells that are killed by effector T cells are assumed to satisfy the Hill equation with Hill coefficient $n=2$ similarly with Eq(\ref{eq13}) the parameter $\beta _{C_{N}}$ represents the maximum killing rate of tumor cell clones. Some tumor clones are assumed to die with the rate $e_{C_{N}}$ due to lack of nutrients and others decays with rate $d_{C_{N}}$ reaching its maximum proliferation level resulting in the action of decomposition.\\ Macrophages participate in essential immune processes in the body and play key roles in human health and disease. The both macrophages phenotypes in Eq. (\ref{eq15}) and Eq. (\ref{eq16}) have a growth rate denoted with $\gamma _{u_{M_{1}}}$ and $\gamma _{u_{M_{2}}}$ support tumor progression by promoting angiogenesis, immune escape, and extracellular matrix remodeling. The interaction of tumor cells and macrophages were also included in Eq. (\ref{eq15}) focusing on the evolution of tumor cells, macrophages with $M_{1}$ like phenotype and macrophages with $M_{2}$ like phenotype. Macrophages also reaches a process were versatile immune cells, shift from one functional state (phenotype) to another, most commonly from pro-inflammatory (M1) to anti-inflammatory (M2) or vice versa, driven by their environment or therapeutic interventions which is also known as polarization of macrophages with rate $\alpha _{m_{2}}$ as modeled in Eq. (\ref{eq16}) with a negative impact \cite{van2016mitochondrial}. \\ In the formulated model we included estrogen receptors (ERs) $E$ in Eq (\ref{eq17}) primarily responsible for normal development, maintaining the extracellular matrix, and regulating elastic recoil in lungs . Estrogen production was modeled as $(1-u)r_{E}F$ were $ r _{E} $ is the estrogen receptors production rate and the parameter $p=1-u$ reduces the effects of $r_{E}$ by $(1-u)r_{E}$ due to aromatase inhibitors in lungs promoting growth of tumor cancer cells. Estrogen are particularly recognized as the pneumocytes and bronchial epithelial cells, which are the primary subtype in lung tissue reaches a degradation level through the ubiquitin-proteasome pathway with rate $\mu$ \cite{maitra2021targeting,durovski2023insights}.The conditions of lipoid pneumonia and lung cancer can coexist and complicate each other through mechanisms like bronchial obstruction, impaired immune function, and the direct spread of cancer cells \cite{maitra2021targeting}. The lipoid pneumonia was modeled in Eq (\ref{eq18}) with a logistic growth having the rate of growth denoted with $r_{F}$. The lung cancer can affect tumors can block airways, leading to poor clearance of secretions, loss of energy with rate $\alpha$ and fats in the lungs, which can cause an increase in lipoid pneumonia.\\ In Eq. (\ref{eq20}), Immature dendritic cells can become mature by capturing free antigens in tissues and antigens on the surface of living tumor cells. In Eq (\ref{eq20}) $\lambda _{D}$ represents the maximum activation rate of dendritic cells, $k_{2}$ is the concentration of antigen that can be recognized on the surface of living tumor cells, $M_{p}$ the half-saturation coefficient of antigens and $d_{D}$ is the natural death rate of mature dendritic cells \cite{xiong2025mathematical, wang2025optimal}. Modeling dendritic cells in lung cancer is an important part of study which are the primary antigen-presenting cells (APCs) responsible for initiating anti-tumor T cell immunity and they are crucial for developing and testing new treatment approaches. Modeling their functions helps in understanding how the immune system naturally recognizes and attempts to fight lung cancer, and why this process often fails in advanced disease states.\\ Lastly in Eq (\ref{eq21}) the activation of effector T cells because of tumor illness through mature dendritic cells by cytotoxic T-lymphocyte antigen . Additional in Eq. (\ref{eq20}) the population of cytotoxic T lymphocytes cells (T) is increased by the dendritic cells which becomes matured through the administered therapeutic particles with $\lambda _{T}$ the maximum activation rate of effector T cells, $M_{T}$ is the half-saturation coefficient of complex cytotoxic T lymphocytes cells and $d_{T}$ is the natural death rate of effector T cells. Effector cells are cleared from the body primarily by phagocytes through a process called apoptosis or programmed cell death once their job of clearing an infection or foreign material is complete but the tumor cell clones is most likely to replace the normal cells their tasks and mechanism resulting in them identifying effector T cells as their threat leading to clearance of effector cells by tumor cells clones with the rate denoted with $\beta _{T}$ also modeled from the Michaelis-Menten kinetics equation.\\ The parameters are also summarized and describe in Table (\ref{T1}) \begin{table}[!h] %\begin{table}[H] \caption{Parameters and values of the ODE model} \label{T1} \small\small \begin{center} \begin{tabular}{l l l l l} \hline \\ Symbol & Description & Units & Value & Source \\ \hline\\ $d_{m_{1}}$ & Natural death rate of $M_{1}$ cells & day$^{-1}$&$0.83-0.924$ & \cite{eftimie2021mathematical}\\ $\alpha _{m_{2}}$& Re-polarization rate of $M_{2}$ macrophages&day$^{-1}$& $0-1$&\cite{eftimie2021mathematical}\\ $\gamma _{C}$ & intrinsic rate of tumor cells $C$&day$^{-1}$&$0\leq \gamma _{C}\leq 1$& \cite{xiong2025mathematical}\\ $\gamma _{C_{N}}$ & intrinsic rate of tumor cell clones & day$^{-1}$&$0\leq \gamma _{C_{N}}\leq 1$& \cite{xiong2025mathematical}\\ $d _{C}$ & Decay rate of sensitive cells&day$^{-1}$&1& \cite{xiong2025mathematical}\\ $n$ & Hill's coefficient&day$^{-1}$&$2$& \cite{xiong2025mathematical}\\ $e_{C_{N}}$ & Death rate of tumor cell clones &day$^{-1}$&$0\leq e _{C_{N}}\leq 1$& \cite{wang2025optimal}\\ $\epsilon$ & Vaccine effectiveness &day$^{-1}$&$0\leq \epsilon \leq 1$&Assumed\\ $\beta _{C}$ & Maximum killing rate & day$^{-1}$&$0 \leq \beta _{C} \leq 1 $ & \cite{xiong2025mathematical}\\ $k_{1}$ & Rate at which antigenic peptides \\&are released by lytic tumor cell &mol.$L^{-1}$.cell$^{-1}$&$3.9664\times 10^{-7}$& \cite{xiong2025mathematical}\\ $F_{Max}$ & Maximum fat growth &$mm^{3}$&$ 0.002711$& \cite{Akman}\\ $u$ &Effectiveness of Aromatase Inhibitors\\&(AI) to limit estrogen production&-&0-1& \cite{Akman}\\ $\mu$ &Estrogen degradation rate& day$^{-1}$ &$1$& \cite{xiong2025mathematical}\\ $\alpha$ & Fat consumption rate &day$^{-1}$mm$^3$&$1.7\times 10^{-6}$& \cite{Akman}\\ $\alpha _{m_{1}}$ & Polarization rate of macrophages \\& with phenotype $M_{1}$ &day$^{-1}$ &$10^{-5}-10^{-2}$&\cite{eftimie2021mathematical}\\ $K_{m}$ & Macrophages carrying capacity & $vol$& $1$ & \cite{eftimie2021mathematical} \\ $d_{P}$ & Degradation rate of antigenic peptide &day$^{-1}$&$1.0082$& \cite{xiong2025mathematical}\\ $\lambda_{D}$ & Maximum growth rate of dendritic cells& day$^{-1}$&$0.1000$& \cite{xiong2025mathematical}\\ $k_{2}$ & Concentration of antigenic peptide\\ & that can be captured on living tumor cells&mol.$L^{-1}$.cells$^{-1}$&$3.9664\times 10^{-8}$& \cite{xiong2025mathematical}\\ $D_{T0}$ & Initial values of dendritic cells &cells&$0$&\cite{xiong2025mathematical}\\ $T_{0}$&Initial value of effector T cells &cells&$0$& \cite{xiong2025mathematical}\\ $d_{T}$ & Death rate of effector T cells & day$^{-1}$&$\geq 1$& Assumed\\ $\beta _{T_{M}}$ & Transmission rates & day$^{-1}$& $\leq 1$ & Assumed\\ $\lambda _{D}$ & maximum activation rate of dendritic cells & day$^{-1}$& $0\leq \lambda _{D}\leq 1$ & \cite{xiong2025mathematical}\\ $\lambda _{T}$ & maximum activation rate of effector T cells & day$^{-1}$& $4$ & \cite{xiong2025mathematical}\\ $\gamma _{u_{M_{1}}}$ & Proliferation rate of macrophages \\ & with phenotype $M_{1}$ & day$^{-1}$& $0.487-0.88$ & \cite{eftimie2021mathematical} \\ $\gamma _{u_{M_{2}}}$ & Proliferation rate of macrophages \\ & with phenotype $M_{2}$ & day$^{-1}$& $0.487-0.88$ & \cite{eftimie2021mathematical} \\ $\lambda _{E}$ & Maximum activation rate of estrogen & day$^{-1}$ & $0-1$ & \cite{van2016mitochondrial}\\ $N_{0}$ & Number of susceptible cells & cells & $\geq 9\times 10^{9}$ & Assumed\\ \hline \end{tabular} \end{center} \end{table} The parameters and initial values of the model with equations from (\ref{eq13}) to (\ref{eq21}) are nonnegative. In fact, the processes ,such as cell activation and division, last from hours to days.The interaction between protein ligands and receptors is of ten completed with in milliseconds to seconds. \subsection{Model Properties} In this section, we demonstrate that the solutions of models from Eq.(\ref{eq13}) to Eq.(\ref{eq21}) exist, are unique and bounded for $t\geq 0.$ Thus, we aim to show that models from Eq.(\ref{eq13}) to Eq.(\ref{eq21}) is mathematically and biologically well-poised. We begin by claiming the following result. \begin{theorem} Model (\ref{eq13})-(\ref{eq21}) with any non-negative initial conditions $(C(0),$ $C_{N}(0),$ $u_{M_{1}}(0),$ $u_{M_{2}}(0), $ $E(0),$ $F(0), $ $P(0), $ $ D_{T}(0), $ $ T_{c}(0))$ $\in\mathbb{R}^6_+$ $:=\{C,C_{N}, u_{M_{1}},u_{M_{2}},E,F,P,D_{T},T_{c},\}:C\geq0,C_{N}\geq 0,u_{M_{1}}\geq0,u_{M_{2}}\geq0,E\geq0,F\geq0\,P\geq0,D_{T}\geq0,T_{c}\geq0 \}$ has a unique solution that is non-negative and uniformly bounded above for all $t\geq0.$\label{theorem1.1} \end{theorem} \begin{proof} We aim to show that if the initial values of System from Eq.(\ref{eq13}) to Eq.(\ref{eq21}) are positive, that is,\\ $C>0$, $C_{N}>0$, $u_{M_{1}}>0$, $u_{M_{2}}>0$, $E(0)>0$, $F(0)>0$, $P>0$, $D_{T}(0)>0$, $T_{c}(0)>0$.\\ From the first equation of system OF equations, we evaluate the derivative at \begin{eqnarray} \ds \left. {\begin{array}{c}\frac{dC}{dt}\end{array} }\right|_{C=0}&=&0\label{eq23} \end{eqnarray} Then also (\ref{eq14}) resulting in the same conditions \begin{eqnarray} \ds \left. {\begin{array}{c}\frac{dC_{N}}{dt}\end{array} }\right|_{C_{N}=0}&=&0\label{eq24} \end{eqnarray} This implies that $C(t)\geq 0$ for all $\forall$ $t\geq \overline{t}$ contradicting the assumption that $C(\overline{t})<0$. Therefore, $C(t)\geq 0$ for all $t>0$. Applying the same reasoning to the other state variables, we obtain \begin{eqnarray} \ds \left. {\begin{array}{c}\frac{dE}{dt}\end{array} }\right|_{E=0}&=&(1-u)r_{E}F\geq 0\label{eq25}\\[5pt] \ds \left. {\begin{array}{c}\frac{dP}{dt}\end{array} }\right|_{P=0}&=&K_{1}\beta _{p}\left( {\begin{array}{c}\frac{T^{2}}{T^{2}+M_{P}}\end{array} }\right)(C+C_{N})\geq 0\label{eq26}\\[5pt] \ds \left. {\begin{array}{c}\frac{dD_{T}}{dt}\end{array} }\right|_{D=0}&=&\lambda _{D}D_{T0} \left( {\begin{array}{c}\frac{P+k_{2}(C+C_{N})}{P+k_{2}(C+C_{N})+M_{P}}\end{array} }\right)\geq 0\label{eq27}\\[5pt] \ds \left. {\begin{array}{c}\frac{dT_{c}}{dt}\end{array} }\right|_{T_{c}=0}&=&\lambda _{T}T_{c0}\left( {\begin{array}{c}\frac{D_{T}}{D_{T}+M_{D}}\end{array} }\right)\geq 0\label{eq28} \end{eqnarray} Thus, we conclude that\\\\ $E(t)\geq 0$, $P(t)\geq 0$, $D_{T}(t)\geq 0$, $T(t)\geq 0$\\\\ This completes the proof. \end{proof} From the Theorem (\ref{theorem1.1}) , when there is no tumor infection in the lungs (that is $\frac{dC}{dt}=0$ and $\frac{dC_{N}}{dt}=0$), the dendritic cells remain in an immature or resting state and perform crucial functions related to maintaining immune tolerance, estrogen levels themselves generally remain within a normal systemic range and self-antigen peptides derived from normal cellular proteins are continuously processed. The macrophages in the absence of tumor continues in maintaining tissue homeostasis through a variety of vital functions, including the routine removal of dead cells, clearance of pathogens and debris, and tissue repair. \newpage %\subsection{Migration Model} %The transportation model or the migration model with PDEs was constructed showing the framework used to analyze and understand the factors that influence the movement of the tumor cell population and body fluids involved in influencing the growth and decomposition of tumor cells. Migration models help to elucidate the intricate mechanisms by which cancer cells move and invade surrounding tissues in the lungs. We first considered the three equations; the first denoted with a variable $C$ was modeled to describe the movements of tumor cells, then the movements of immune cells denoted with a variable $n$ which are activated in response of infected normal cells and then the movement of microscopic particles from chemotherapy drugs and any other form of treatment used denoted with $T_{D}$ providing a combined directed treatment response of tumor in lungs \cite{palucka2012cancer}. \\ %Cancer cells and their associated components take advantage of the body's many fluid systems and their properties ( eg temperature, viscosity etc) as a means of dissemination throughout the body and colonization of distant organs \cite{tarin2011cell}. The motion of the body fluids especial those in lungs are considered in the model and denoted with $v$ representing the fluids velocity with applications of the continuity equation in PDEs form with the mass flow rate of fluids in lungs but taken as constant within a system, which is based on the principle of conservation of mass. Cancerous tissues often have higher metabolic activity and increased blood flow compared to healthy tissue, resulting in elevated surface temperatures $T$. Temperature modeling helps interpret these thermal activities to accurately detect and locate tumors, particularly in the applications of lung cancer \cite{mital2007thermal}. The blood flow also follows the Lagrangian approach a method which focuses on individual particle through the flow field and the fluid properties are given for a fixed particle for varying time \cite{shadden2015lagrangian}. Lagrangian is used in experiments where a trace is thrown into moving fluid bodies for to determine currents for example when a pharmaceutical compound or a drug is being observed its movements from the point of release to the targeted area the blood where most of tumor cells are located \cite{shadden2015lagrangian,tarin2011cell}.\\ %Applying starling principle the fluid movement across a semi-permeable blood vessel such as a capillary or small venule located in lungs (where gas exchange occurs) is determined by the hydrostatic pressures and colloid osmotic pressures (oncotic pressure) on either side of a semipermeable barrier that sieves the filtrate, retarding larger molecules such as proteins from leaving the blood stream. As all blood vessels allow a degree of protein leak , true equilibrium across the membrane cannot occur and there is a continuous flow of water with small solutes \cite{curry1980fiber}. The molecular sieving properties of the capillary wall reside in a recently discovered endocapillary layer rather than in the dimensions of pores through or between the endothelial cells. From these applications then the transportation model satisfies the following model in Eq. (\ref{eq32}) to Eq. (\ref{eq36}) %\begin{eqnarray} %\ds \underbrace{\frac{\partial C}{\partial t}+\frac{1}{r}\frac{\partial}{\partial r}(r^{2}uC)-\delta _{C}\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial C}{\partial r}\end{array} }\right)}_\textrm{Tumor cells concentration velocity equation}&=&\ds \underbrace{\gamma _{C}C \left( {\begin{array}{c} 1-\frac{c}{K_{C}}\end{array} } \right)}_\textrm{Logistic growth term}-\underbrace{\frac{4k_{0}}{3r^{2}}\frac{\partial }{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial (T^{4})}{\partial r}\end{array} }\right)}_\textrm{$ \left. \begin{tabular}{ll} %The skin temperature or\\ %radiative heat transfer % \end{tabular} \right. $ }\cr &&+\underbrace{\frac{c_{g}\rho _{b}}{4\pi K}\left[ {\begin{array}{c}\frac{(1-\frac{\rho _{p}}{\rho _{T}})^{3}}{(\frac{\rho _{p}}{\rho _{T}})^{2}}\end{array} } \right]\frac{M_{c}av}{F}C}_\textrm{$ \left. \begin{tabular}{ll} %Volume of toxic gases with \\ %radon level from cigarette %\end{tabular} \right. $}-\underbrace{\beta _{C}\frac{nC}{n+g}}_\textrm{$ \left. \begin{tabular}{ll} %Clearance of tumor by \\ %immune cells %\end{tabular} \right. $}\label{eq32}\\[5pt] %\ds \underbrace{\frac{\partial n}{\partial t}+\frac{1}{r}\frac{\partial}{\partial r}(r^{2}un)-\delta _{n}\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial n}{\partial r}\end{array} }\right)}_\textrm{Immune cells response velocity equation}&=&\ds %\underbrace{\lambda _{n}n(1-n)}_\textrm{Proliferation of immune cells}-\underbrace{\beta _{n}\frac{nC}{C+g}}_\textrm{$ \left. \begin{tabular}{ll} %Clearance of immune \\ %cells by tumor cells %\end{tabular} \right. $ }\cr && %+\underbrace{\frac{4k_{0}}{3r^{2}}\frac{\partial }{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial (T^{4})}{\partial r}\end{array} }\right)}_\textrm{$ \left. \begin{tabular}{ll} %The skin temperature or\\ %radiative heat transfer %\end{tabular} \right. $ }+\underbrace{\chi C_{n}}_\textrm{$ \left. \begin{tabular}{ll} %sensitivity of immune \\ %cells to chemokines %\end{tabular} \right. $ }\label{eq33}\\[5pt] %\ds \underbrace{\frac{\partial C_{T}}{\partial t}+\frac{1}{r}\frac{\partial}{\partial r}(r^{2}uC_{T})-\delta _{C_{T}}\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial C_{T}}{\partial r}\end{array} }\right)}_\textrm{Therapeutic particles velocity equation}&=&\ds \underbrace{\eta _{0}\left[{\begin{array}{c}\frac{2}{r}\frac{\partial n}{\partial r}+\frac{\partial ^{2} n}{\partial r^{2}}\end{array} }\right]}_\textrm{$ \left. \begin{tabular}{ll} %Direct boost of immune cells \\ %by therapeutic particles %\end{tabular} \right. $}-\underbrace{\alpha _{C_{T}}C_{T}}_\textrm{$ \left. \begin{tabular}{ll} %Wash out %\end{tabular} \right. $ }\cr &&+\underbrace{\beta _{C_{T}}\frac{nC}{g+n}\left[ {\begin{array}{c}\frac{1}{1+Q/K_{TQ}}\end{array} }\right]}_\textrm{$ \left. \begin{tabular}{ll} %Activation of immune cells\\ %and inhibition by PD-1 % \end{tabular} \right. $ }-\underbrace{d_{C_{T}}C_{T}}_\textrm{ $ \left. \begin{tabular}{ll} %metabolic\\ %degradation %\end{tabular} \right. $ }\label{eq34}\\[5pt] %\ds \underbrace{\frac{\partial v}{\partial t}+\frac{1}{r}\frac{\partial}{\partial r}(r^{2}uv)-\delta _{v}\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial v}{\partial r}\end{array} }\right)}_\textrm{Blood velocity equation}&=& \ds \underbrace{-\frac{1}{\rho _{f}}\frac{\partial P}{\partial r}}_\textrm{ $\left. \begin{tabular}{ll} %Blood pressure gradient\\ %force per unit mass %\end{tabular} \right. $ }+\underbrace{\frac{c_{b}\rho _{b}}{k_{1}}\frac{\partial T}{\partial t}}_\textrm{Heat conduction}\cr && % -\underbrace{\frac{\rho _{b}c_{b}w_{b}}{k_{1}}(T_{a}-T_{b})}_\textrm{Blood perfusion term}+\underbrace{L_{p}S(\Delta P_{i}-\sigma \Delta \pi _{i})}_\textrm{$\left. \begin{tabular}{ll} %Pressure on semi-\\ %permeable barrier %\end{tabular} \right. $}\label{eq35}\\[5pt] %\ds \underbrace{\frac{\partial T}{\partial t}+\frac{1}{r}\frac{\partial}{\partial r}(r^{2}uT)-\delta _{T}\frac{1}{r^{2}}\frac{\partial}{\partial r}\left( {\begin{array}{c}r^{2}\frac{\partial T}{\partial r}\end{array} }\right)}_\textrm{Temperature for fluids in motion equation}&=&\ds \underbrace{\frac{v}{c_{p}}\left( {\begin{array}{c}\frac{\partial v}{\partial t}\end{array} }\right)^{2}}_\textrm{Energy transfer process}+\underbrace{\rho _{b}c_{b}w_{b}(T_{a}-T_{b})}_\textrm{Blood perfusion term}\cr &&+\underbrace{\tau \left[ {\begin{array}{c}D_{B}\frac{\partial C}{\partial t}\frac{\partial T}{\partial t}+\frac{D_{T}}{T_{\infty}}\left( {\begin{array}{c}\frac{\partial T}{\partial t}\end{array} }\right)^{2}\end{array} }\right]}_\textrm{ $\left. \begin{tabular}{ll} %Temperature variation due\\ %to tumor concentration %\end{tabular} \right. $ }\label{eq36} %\end{eqnarray} %In both of the models in Eq.(\ref{eq32}) to Eq. (\ref{eq36}) the equation on the left hand side given represents the equation for velocity of particles. The right hand side of Eq. (\ref{eq32}) , tumor cell velocity follows a logistic growth as in Eq. (\ref{eq13}). Then also the next equation on the right hand side of Eq. (\ref{eq32}) which is modeled to show the positive effectiveness of radiative heat flux described by Stefan Boltzmann law with a $k_{0}$ as the Boltzmann's constants which is thermal physiology brunch that deals %with normal functions of living organism stating that the skin temperature is a resultant balance of the heat transport within the tissue and its transport to the environment, the equation is also applicable in the response of immune cells given in Eq. (\ref{eq33}) with a slightly negative impact \cite{jasinski2010influence}.\\ %In the lungs the sensitivity of immune cells to chemokines $\chi$ with sensitivity of concentration of chemokines $C_{n}$ is dynamically regulated by the local microenvironment, particularly during inflammation or infection, which generally increases their sensitivity and directs their movement. The flow of smoke through the tobacco rod and filter is generally modeled as in Eq.(\ref{eq32}) applying the Darcy’s Law for a steady or unsteady laminar flow of an incompressible fluid (or sometimes compressible) fluid through a porous medium \cite{mcadam2016influence} with parameters $K$ denoting the empirical permeability factor, and $c_{g}$ denoting circumference of the cigarette. The volumetric flow rate of smoke ($F$) through the cigarette is related to the pressure drop $P$ with the total mass flow rate as a product of the fluid density ($\rho$ ), cross-sectional area ($A$) also given as $A=\frac{F}{P}$, and velocity $v$.\\ %In Eq. (\ref{eq34}) the therapeutic particles especially when combined with nanoparticles activates the immune system with a positive reaction rate $\eta _{0}$ to treatment process but generally reduces their overall circulation time and leads to their eventual removal from the bloodstream with rate $\beta _{n}$ \cite{shadden2015lagrangian}. In living tissue, blood circulation plays an important role with a duel purposes, it allows nutrient and gas exchanges which are crucial for the organic function while serving as a natural temperature regulator \cite{ndreko2025modeling}. In order to incorporate this processes the bio-heat equation in Eq. (\ref{eq35}) and Eq. (\ref{eq36}) was introduced in the model by incorporating the perfusion term with $\rho _{b}$ as blood density, $c_{b}$ blood heat capacitance and $k_{1}$ tissue thermal conductivity \cite{ndreko2025modeling}. The classic case of the heat conduction is derived from the original heat equation in a homogeneous medium (a lung branch level in this case) with $T$ defining the temperature and as conduction only goes on x-axis.\\ %The application of Starling principle was considered in Eq. (\ref{eq35}) which holds that fluid movement across a semi-permeable blood vessel such as a capillary or %small venule is determined by the hydrostatic pressures $P_{i}$ and colloid osmotic pressures (oncotic pressure) $P_{c}$ on either side of a semipermeable barrier that sieves the filtrate, retarding larger molecules such as proteins from leaving the blood stream. The therapeutic particles modeled in Eq. (\ref{eq35}) is assumed to be washed out or cleared from the human body with rate $a_{C_{T}}$,which generally signifies the end of their therapeutic action, followed by metabolic degradation with rate $d_{C_{T}}$ and excretion. The primary consequence of this elimination is the cessation of the drug's therapeutic effect, necessitating further doses to maintain treatment. \\ %The molar flux shown at the right hand side of Eq. (\ref{eq36}) associates species in the blood vessels transfer by diffusion determined by an expression from Fourier’s law for the conditions which is termed to Fick's law with the equation $D_{B}\frac{\partial C}{\partial t}$ where $D_{B}$ is a property of the binary mixture known as the binary diffusion coefficient. The parameters are also summarized in Table (\ref{T2}) along with their units and values. The migration model will help understanding the complex internal signals, cells migrations, treatment resistance by explaining how certain cell clusters survive in the bloodstream and resist immune surveillance and therapies, contributing to more successful metastatic events. However the primary goal of studying migration is to find ways to stop it through crucial tools for developing and testing new cancer treatments like ; drug testing; target identification and ranking cancer origin \cite{welf2011signaling, sengupta2021principles}. %\begin{table}[!h] %\begin{table}[H] %\caption{Parameters and values for the migration model} %\label{T2} %\small\small %\begin{center} %\begin{tabular}{l l l l l} %\hline \\ %Symbol & Description & Units & Value & Source \\ %\hline\\ %$k_{0}$ & The Boltzmann constant & $J/K$&$1.380649 \times 10^{-23}$& Assumed\\ %$Q$ & Inhibition by the complex PD-1\\ & an immune checkpoint protein & percentage ($\%$) &$0\leq Q\leq 100$& Assumed\\ %$c_{b}$ & Blood heat capacitance &$mL/mmHg$&$2\leq c_{b}\leq 4$&\cite{ndreko2025modeling}\\ %$\rho_{b}$& Blood density &$g/mL$&$1.05\leq \rho _{b} \leq 1.06$& Assumed\\ %$L_{p}$ & Hydraulic conductivity of the membrane & $cm.mm.Hg^{-1}.s^{-1}$&$5.0\times 10^{-8}$& Assumed\\ %$S$ & The surface area for filtration & units$^{2}$& $\leq 1$ & Assumed\\ %$\Delta P_{i}$ & Difference between capillary hydrostatic \\& and interstitial hydrostatic pressure & $mmHg$ & $\geq 1$ & Assumed\\ %$\Delta \pi _{i}$ & Difference between plasma protein oncotic \\& and subglycocalyx oncotic pressure & $mmHg$ & $\geq 1$ & Assumed\\ %$\rho_{p}$ & Packing density of tobacco & $g/cm^{3}$ & $0.287$ & \cite{mcadam2016influence}\\ %$\rho_{T}$ & Density of tobacco shreds & $g/cm^{3}$ & $0.25$ to $0.287$ & \cite{mcadam2016influence}\\ %$\sigma$ & Staverman's reflection coefficient & $mL.mm.Hg^{-1}.s^{-1}$ & $\geq 1$ & Assumed \\ %$T_{\infty}$ & Temperature at point & $K$ & $\geq 15$ & Assumed \\ %$D_{B}$ & Diffusion coefficient & $cm^{2}s^{-1}$ & $\leq 10^{-5}$ & \cite{ndreko2025modeling}\\ %$w_{b}$ & Blood perfusion rate & $mL(min.g)^{-1}$ & $>5.0$ & \cite{ndreko2025modeling}\\ %$K$ & Empirical permeability factor & $m^{2}$ & $9.86923\times 10^{-13}$ & \cite{mcadam2016influence}\\ %$F$ & Flow rate of smoke & $cm^{3}s^{-1}$ & $0-1$ & \cite{mcadam2016influence}\\ %$M_{c}$ & Molar mass of gas with radon level & $g/mol$ & $222.01758$ & \cite{mcadam2016influence}\\ %\hline %\end{tabular} %\end{center} %\end{table} \subsection{Drug Administration} To model the dynamic drug administration process, we let the administration method be oral where the therapeutic particles are taken by mouth and intravenous which is common for chemotherapy and immunotherapy given into a vein via needle or catheter over time, often in a clinic or hospital \cite{cummings2008administration, holliday2008administration}. Let $j$ (the $(j-1)$'th day) represent the potent time points of drug intake, and $N$ denotes the total number of therapy days \cite{wang2025optimal}. According to Wang effective dose of the drug (cisplatin, carboplatin, gemcitabine, paclitaxel, and docetaxel) $D (t)$ in the human body at time t can be formulated as \begin{eqnarray} \ds \frac{dD_{j}}{dt} &=& (\varphi_{0}-\alpha _{r})\sum_{j=1}^{N}D _{j}(t) \label{eq37} \end{eqnarray} where \begin{eqnarray} \ds D(t) &=& \left \{ {\begin{array}{c c} D _{j,0}e^{(\varphi_{0}-\alpha _{r})\Delta t_{j}} &, t_{j} \leq t \leq t_{j}+1,\\\\ 0 &, \textrm{otherwise}~~~~~~~\\ \end{array} }\right. \label{eq38} \end{eqnarray} Here, $\beta _{i,0}$ represents the administered dose at $t=t_{i}$ arranging from $0$ to the maximum tolerated dose (MTD) M, $a_{r}=\alpha _{0}+\mu _{\beta}$ represents the drug resistance rates with $\alpha _{0}$ resistant rate due to body temperatures and $\mu _{\beta}$ resistant rate due to decay \cite{xiong2025mathematical, wang2025optimal} and $\varphi_{0}$ is the drug administration rate which must be less than the resistance rate ($\varphi\leq a_{0}$) for a decrease in concentration. The drug resistance rate $a_{r}$ shown in Eq. (\ref{eq37}) is positively influenced by body temperature, density $\rho _{b}$ of blood due to increased blood viscosity $\mu$ , changes in the diameter of blood pathways in veins causing obstruction of blood flow due to tumor volume $V$ in the lungs, and also energy lost $h_{f}$ due to the consumption of fat by the tumor. Using Buckingham's pi theorem simplifying complex physical problems by reducing the number of variables into dimensionless groups revealing underlying physical relationships to analyze the conditions of drug resistance $a_{r}$ we get Eq. (\ref{eq39}) \begin{eqnarray} \ds a_{r} &=& \frac{M}{\Delta t}R_{e}f( \Pi _{s}) \label{eq39} \end{eqnarray} where $f(\Pi _{s})=f \left ( {\begin{array}{c c c} R_{e} & , \frac{V}{D^{3}}&,\frac{h_{f}}{\rho U D^{4}} \end{array} }\right)$ a function of Reynold's number $R_{e}$ given by $R_{e}=\frac{\mu}{\rho U D}$ and of Teichholz formula to account for the changing width-to-length ratio of veins in lungs. $M$ represents the mass of drug microscopic particles and $\Delta t$ represents the change in time at each point of drug resistance together Eq. (\ref{eq39}) is the rate of change of mass describing the drug losing mass over time. The Teichholz formula is widely used in echocardiography to estimate the Left Ventricular (LV) Volume and Ejection Fraction (LVEF) based on linear measurements from a one-dimensional $2D$-derived M-mode. In the lungs, echocardiography acts as a non-invasive ultrasound-based imaging test used to assess how lung diseases affect the right side of the heart, particularly in diagnosing pulmonary hypertension, blood clots (pulmonary embolism), or assessing the impact of lung cancer disease on heart functionality which results in reduced effectiveness of drugs on the disease. \section{Parameter Estimations and Data Fittings} Our mathematical model aims to quantify the ability of an intervention, drug, or action to produce a desired, beneficial result under ideal and controlled circumstances of the lung cancer drugs on tumor size. It serves as a framework for combining the effects of lung cancer therapeutic particles . In the following, we estimate the specific initial conditions and parameter values for modeled equations in Eq. (\ref{eq13}) to Eq. (\ref{eq21}). \subsection{Parameter Estimations} Parameters were categorized into public parameters, which remains constant across patients and are derived from literature and fitting data, and individualized parameters, which vary between patients due to the heterogeneous nature of lung cancer response to therapy. To determine the number of cells based on tumor size, we consider the tumor to be spherical and the diameter of each cell is $20\mu m$ \cite{xiong2025mathematical}. In Liu et al. (2018), the tumor size unit is $mm^{3}$ \cite{lu2018haart}. Therefore, the relationship between tumor size and the number of tumor cells is determined by the following equation \begin{eqnarray} N_{0}(L)=\frac{R_{b}L}{\frac{4\pi}{3}\left ( {\begin{array}{c} \frac{0.02}{2}\end{array} }\right )^{3}\label{eq40} } \end{eqnarray} Where $L$ is the tumor size and $R_{b}=0.7405$ is a spherical packing coefficient \cite{xiong2025mathematical}. The parameters of Models in Eq.(\ref{eq13}) to Eq.(\ref{eq21}) are estimated as follows. \begin{itemize} \item [(i)] {\bf Semi-Saturation Coefficients}\\ The semi-saturation coefficients can be estimated by the steady-state approximation method \cite{lai2018modeling}. An expression of the form $\frac{A}{A+K_{A}}B$ indicates that $B$ is activated by $A$, the half-saturation coefficient $K_{A}$ is considered as the approximate steady-state concentration of species $A$ \cite{xiong2025mathematical}. Applying the method to estimate the half saturation coefficients, we assume that the steady state of System for Eq.(\ref{eq13}) to Eq.(\ref{eq21}) without treatments is ($E_{C}$, $E_{C_{N}}$, $E_{u_{M_{1}}}$, $E_{u_{M_{2}}}$, $E_{E}$,$E_{F}$,$E_{P}$,$E_{D_{T}}$,$E_{T_{c}}$). In a condition that the systems of equations in Eq.(\ref{eq13}) to Eq.(\ref{eq21}) are in a steady state the macrophages turn to be constant meaning they are persistently present and active in a tissue rather than arriving, performing a temporary cleaning function, then we assume that the two macrophages are equal $u_{M_{1}}=u_{M_{2}}$ and also the tumor cells becomes constant thus the tumor's size, volume, or growth rate is treated as unchanging over a specific period making the two tumors equal $C=C_{N}$. Then for the Semi saturation coefficients is reduced from the condition of stead state equating the Eq.(\ref{eq15}) to zero. For the form $\frac{C+C_{N}}{C+C_{N}+M_{u_{M_{1}}}}$, based on a steady state approximation we have $M_{u_{M_{1}}}=E_{C}+E_{C_{N}}$, then \begin{eqnarray} \ds \gamma _{u_{M_{1}}}E_{u_{M_{1}}}\left( {\begin{array}{c} 1-\frac{E_{u_{M_{1}}}}{K_m}\end{array} } \right)-\frac{1}{2}\alpha _{m_{1}}E_{u_{M_{1}}}-d_{m_{1}}E_{u_{M_{1}}}&=&0\label{eq42} \end{eqnarray} Therefore, we have a steady state condition for macrophages neglecting the natural death, since it appears to be a silent process occurring without triggering inflammation \begin{eqnarray} E_{u_{M_{1}}}\leq K_{m}\left( {\begin{array}{c}1-\frac{\alpha _{m_{1}}}{2\gamma _{u_{M_{1}}}}\end{array} } \right)\label{eq43} \end{eqnarray} For the steady state of tumor cells the natural decay (death) rate of tumor cells becomes equal to their proliferation (birth) rate, resulting in no net growth. Instead of decreasing or being eliminated, natural death and cell division continue at a balanced, constant rate, maintaining a stable tumor volume. Simplifying Eq.(\ref{eq13}) in a steady state form using the described conditions gives the following equation: with $\omega _{0}=\beta _{C}+\beta _{D}$, \begin{eqnarray} \ds E_{C}&\geq&K_{C}\left( {\begin{array}{c}1-\frac{\omega _{0}}{2\gamma _{C}}\end{array} } \right)\label{eq44} \end{eqnarray} Similarly, for the form $\frac{P+k_{2}(C+C_{N})}{P+k_{2}(C+C_{N})+M_{P}}$ the steady-state approximation method is used again to obtain $M_{P}=E_{P}+2k_{2}E_{C}$. Then, when fat availability is constant during a tumor infection, the tumor typically hijacks normal fat metabolism to accelerate its own growth, promoting metastasis, and involuntary weight loss in the host; this reduces the equation of modeled fats in Eq.(\ref{eq18}) with steady state condition to become $E_{F}=F_{Max}$ \cite{martin2022role}. In summary, the estimated half-saturation coefficients are \begin{eqnarray} \ds &&M_{u_{M1}}\geq2K_{C}\left( {\begin{array}{c}1-\frac{\omega _{0}}{2\gamma _{C}}\end{array} } \right),~~~~~M_{P}\geq K_{C}\left( {\begin{array}{c}1-\frac{\omega _{0}}{2\gamma _{C}}\end{array} } \right)\left( {\begin{array}{c}\frac{k_{1}\beta _{P}}{d_{P}}+2k_{2}\end{array} } \right) \cr && M_{D}=\frac{\lambda _{D}D_{T0}}{2d_{D}}, ~~~~~M_{C_{N}}\geq K_{c}\left( {\begin{array}{c}1-\frac{\omega _{0}}{2\gamma _{C}}\end{array} } \right), ~~~~~~M_{T_{C}}=\frac{\lambda _{T_{c}}T_{0}}{\beta _{T_{0}}-2d_{T_{c}}}\cr && M_{E}=\frac{\lambda _{E}E_{0}}{4\mu}+\frac{r_{E}}{\mu}(1-u)F_{Max} ~~~or~~M_{E}=2K_{m}\left( {\begin{array}{c}1-\frac{\alpha _{m_{1}}}{2\gamma _{u_{M_{1}}}}\end{array} } \right)\label{eq45} \end{eqnarray} \item [(ii)] {\bf Maximum tumor cell growth rate}\\ In Xiong, Z at al (2025) \cite{xiong2025mathematical} cancer cells are suitable for experimental animal models to study human lung cancer for its high tumorigenic potential and aggressive nature. Lots of experiments have shown that the doubling time of cancer cells is $13.6 \pm 1.5$ hours in an ideal vitro condition \cite{simoes2015metabolic}. The exponential model is the natural description of the early stages of cancer growth, where each cancer cell was divided into two daughter cells in the affected area at the constant rate $\gamma _{C}$ from the origin and the rate of tumor clones $\gamma _{C_{N}}$. The estimation of the maximum growth rate parameter comes from the power law model that states that growth is proportional to the volume of tumor cells $C$ as stated in Eq.(\ref{eq13}) but in summary the origin of the equation is as follows, \begin{eqnarray} \ds \frac{dC}{dt}&=&\gamma _{C}C^{\alpha}\label{eq46} \end{eqnarray} where $\gamma _{0}$ is the maximum growth rate and the unit of time $t$ is day using $\alpha=1$. Most lung cancer cases are diagnosed later, making it difficult to measure the speed with which they progress. This is because many people don’t experience symptoms until the cancer has grown and spread. People with non small cell lung cancer are frequently diagnosed at stage $3$ or $4$, but tumors smaller than $4cm$ are classified as stage $1$, those between $4cm$ to $7cm$ as stage $2$, and those larger than $7cm$ as stage $3$ lung cancer \cite{jia2020study}. At these stages, the cancer has typically spread to nearby lymph nodes, nearby lung tissues, other chest structures, and other parts of the body, known as metastasis. According to Jia B ($2020$) \cite{jia2020study} tumor volumes for staging can range from $\leq 2.80cm^{3}$ at stage $1$ up to $\geq 55cm^{3}$ at stage $5$ . \begin{eqnarray} \ds \gamma_{C}&=&\frac{\ln \left( {\begin{array}{c}\frac{C_{b}}{C_{a}}\end{array} } \right)}{\Delta t}=\frac{\ln \left( {\begin{array}{c}\frac{55}{2.80}\end{array} } \right)}{261-70}\geq 0.0156 \label{eq47} \end{eqnarray} Therefore the maximum growth rate of tumor cells $\gamma _{C}\geq0.0156$ day$^{-1}$, but for tumor cell clones applying the same concept in Eq.(\ref{eq47}) from Eq.(\ref{eq46}) the maximum growth rate $\gamma _{C_{N}}\geq0.2221$ day$^{-1}$ which shows a higher progression rate of clonal tumor cells, clinically significant progression is marked by tumors larger than $3 cm$ (Stage $2$) or high-risk volumes exceeding $17,010 mm^{3}$ \cite{ding2025characterising}. \item[(iii)] {\bf Carrying capacity of tumor cells ($K_{C}$ and $K_{C_{N}}$)}\\ The carrying capacity $K$ of tumor cells in the lungs represents the maximum tumor size or volume that lung tissue can support, influenced by the supply of nutrients, vascularization, and the immune response. Two patients with identical tumor volumes can have vastly different growth dynamics according to their individual carrying capacities and therefore different responses to therapy \cite{harshe2023predicting}. The actual carrying capacity is not fixed and varies significantly between patients, determined by factors such as the rate of angiogenesis (blood vessel formation) and the host immune response. Mathematical models, such as the Gompertz model that uses $3.1\times 10^{12}$ as the standard theoretical upper limit often used for many solid tumors \cite{talkington2015estimating}. The carrying capacity of tumor clones $K_{C_{N}}$ are specific, often smaller, niche constraints faced by different sub-clonal populations within the tumor, which can differ based on their phenotype, location, and ability to compete for resources \cite{gomez2020heterogeneity}. \item[(iv)] {\bf Naive mature dendritic cells $D_{T0}$} \\ Dendritic cells are usually not abundant at tumor sites, but increased densities of populations of dendritic cells have been associated with better clinical outcome, suggesting that these cells can participate in controlling cancer progression. Lung cancers have been found to include four different subsets of dendritic cells classified as dendritic cell subsets and one plasmacytoid dendritic cell subset \cite{dempsey2015distinct}. Dendritic cells (DCs) generally take $1$ to $2$ days to fully mature after receiving an activation signal. In vitro, the process of generating mature DCs from monocytes typically takes about $7$ days \cite{nair2012isolation}. Naïve dendritic cells originate from progenitor cells in the bone marrow and migrate to almost all tissues throughout the body \cite{xiong2025mathematical}. According to Eq.(\ref{eq20}) in the model, the stabilized replenishment rate of dendritic cells applying the steady state condition the natural death is described by $\frac{\lambda _{D}D_{T0}}{2}$. In a steady state, the replenishment rate of approximately $4000$ cells per hour is equals to the death rate. Hence we have, \begin{eqnarray} \ds \frac{\lambda _{D}D_{T0}}{2}=4000\times 24 \label{eq48} \end{eqnarray} Then, we solve Eq.(\ref{eq48}) for the naïve dendritic cells to get $D_{T0}= 9.6\times 10^{4}$ cells. \item[(v)] {\bf Naive effector T cells $T_{0}$}\\ In vitro condition the total number of thymocytes or the immature T-cells precursors that originate in the bone marrow and migrate to the thymus (the upper chest between the lungs) to undergo maturation in $5.5$ week old mice is approximately $2\times 10^{6}$ cells \cite{xiong2025mathematical}. Some studies indicates that approximately $0.7$ to $1\%$ of the total thymocytes population is exported to the peripheral immune system as new T cells daily \cite{berzins1998role}. According to the fourth equation in System Eq.(\ref{eq21}), the stabilized replenishment rate of T cells is described by $\frac{\lambda _{T}T_{0}}{4}$. Hence, we set the daily migration rate \begin{eqnarray} \ds \frac{\lambda _{T}T_{0}}{4}&=&\frac{T_{c}}{2}(\beta _{T_{c}}+2d_{T_{c}})\leq(2\times 10^{6})\times (8.0\times 10^{-3})\label{eq49} \end{eqnarray} Then, solving for $T_{0}$, we obtain $T_{0}=1.6\times 10^{4}$ cells. However using Eq. (\ref{eq49}) killing rate of T cells by tumor clones in vitro studies with high efficiency killing achieves $99.9\%$ clarence per day, then \begin{eqnarray} \ds 2&=&2d_{T_{c}}+0.999\label{eq50} \end{eqnarray} Solving Eq.(\ref{eq50}) by equating $\frac{T_{c}}{2}(\beta _{T_{c}}+2d_{T_{c}})$ with the condition that at steady state meaning not yet encountered their specific antigen in the periphery giving $T_{c}(0)=T_{0}$ then $T_{c}=T_{0}$ to get the death rate of effector T cells $d_{T_{c}}=0.5$. Naive T cells are mature, antigen inexperienced T lymphocytes that circulate between blood and lymphoid tissues, they occurs when T Cell Receptors binds to a specific antigens. \item[(vi)] {\bf Estrogen receptors $E_{0}$} \\ The number of estrogen receptors (ER) is highly variable, depending on the cell line, culture conditions (serum, hormone depletion), and passage number. The most commonly studied in vitro model for estrogen receptors reports indicate a range of $1,0\times 10^{4}$ to $5.0\times 10^{4}$ receptors per cell \cite{durovski2023insights}. Patients with estrogen receptors (ER) negative lung cancer have a significantly higher mortality rate. Estrogen receptor beta (ER$\beta$) is the predominant form in lung cancer (detected in up to $60\%$ of cases), whereas ER$\alpha$ is less commonly found in the lung, though it may be present. Applying the steady state condition in Eq.(\ref{eq17}) then we get the following, \begin{eqnarray} \ds \frac{\lambda _{E}E_{0}}{4}&=&\mu E-(1-u)r_{E}F_{Max}\leq (5.0\times 10^{4})\times (9.0\times 10^{-3})\label{eq51} \end{eqnarray} Solving Eq.(\ref{eq51}) to get $E_{0}=4.5\times 10^{2}$ cells applying the steady condition with the maximum law fats consumption approximately equal to zero $F_{Max}\approx 0$. When fat availability approaches zero, tumor cells face severe metabolic stress, leading to a state of starvation, reduced proliferation, and, in many cases, cell death \cite{hsu2017estrogen}. \item[(vii)] {\bf Macrophages with phenotypes $M_{1}$ and $M_{2}$}\\ In vitro condition an experiment was done in \cite{patel2017fate}, which focused on the adoptive transfer of human mononuclear phagocytes into mice, showed that classical ($M_{1}$) macrophages circulate for a mean of $1.01$ days and ($M2$) macrophages have longer mean lifespans of $7.41$ days. The authors in \cite{italiani2014monocytes} also stated that murine $M_{1}$ like macrophages (involved in phagocytosis) have a half-life between $18-20hrs$, while the murine $M2$ like macrophages (involved in tissues repair) have a half-life between $5-7$ days. We assumed that in Eq(\ref{eq13}) the base line value of $r_{m_{2}}$ is $0.1$ such that $r_{m_{2}} \epsilon (10^{-2}-10^{0})$. To avoid numerical problems caused by a stiff system, we rescale both tumor cells (tumor cell and tumor cell clones ) by $K_{C}$ and $K_{C_{N}}$ along with macrophages by $K_{m}$ leading unitary carrying capacities ($\overline{K_{C}}=1, ~~\overline{K_{C_{N}}}=1, ~~\overline{K_{m}}=1$), and to rescale the following parameters: \begin{eqnarray} \ds && \overline{d_{C}}=d_{C}\times K_{m},~~~~ \overline{d_{C_{N}}}=d_{C_{N}}\times K_{m}\cr && \overline{e_{C_{N}}}=e_{C_{N}}\times K_{M}, ~~~~\overline{K_{C}}=\frac{K_{C}^{*}}{K_{C}}, ~~~~~\overline{K_{C_{N}}}=\frac{K_{C_{N}}^{*}}{K_{C_{N}}}\label{eq0} \end{eqnarray} The values in Eq.(\ref{eq0}) are used for the numerical simulations performed throughout this study. Also in the study noted in \cite{patel2017fate} the mixed phenotype are said to have a small potion which nearly approximately to zero and can also be neglected for a simple model. Thus, we assume that when the tumor is introduced into the system, it elicits a pro-inflammatory immune response characterized only by the presence of a non-zero $u_{M_{1}}$ population. The numerical solution is propagated in time using a classical fourth order Runge-Kutta method with the conditions which are numerical simulated in \cite{eftimie2021mathematical}, then \begin{eqnarray} \ds && u_{M_{1}}(0)=\frac{0.006}{K_{m}}=\frac{0.006}{6.72}= 0.000893 \cr && u_{M_{2}}(0)= 0.0 \end{eqnarray} The parameter value of death rate is estimate to be $d_{C}=0.87$ per day assuming that they have a correspondence of $M_{1}$ like macrophages have a half-life of $\approx 0.8$ days, same with $M_{1}$ like macrophages have a half life of $\approx 5.14$, corresponding to a death rate of $0.09$ per day. \end{itemize} \subsection{Data Fitting} When lung cancer grows without control or treatment, as graphically shown in Fig (\ref{Fig:3 Graphical prisanation of tumor growth without control}), it continues to multiply, forming tumors that destroy healthy lung tissue, block airways, and spread to other parts of the body. This progressive appears to be linear as shown in Fig (\ref{Fig:3 Graphical prisanation of tumor growth without control}) if left unchecked, leading to severe growth, worsening symptoms, and, if left entirely untreated, is generally fatal \cite{ozkose2021fractional}. If the intrinsic rate $\lambda _{C}$ becomes more than $3$, then the radiance of the gasses inhaled by patience will be more deadly than smoking cigarettes, which can block the airways in the lungs in an instance or a short period of time. Mostly, this type of growth occurs in a short period of time ($4$ months) if not treated to stage 4 of lung cancer, defined by most as the advanced form of the disease, and is often the most difficult diagnosis for patients and their families to process\cite{proctor2012history}. \\ At this stage, cancer spreads beyond the lungs to distant organs such as the brain, liver, bones, or adrenal glands, making curative treatment extremely challenging \cite{dwivedi2020perspective}. When a diagnosis is referred to as terminal lung cancer or stage $4$, the most immediate and difficult question often centers on life expectancy for those who choose not to seek medical intervention. When the immune system interacts with the tumor, they fight to reduce the growth rate of the tumor. Treatments can also stimulate immune cells for a quicker response against tumor infections resulting in controlled tumor growth. In the early stages, the immune system recognizes tumor specific antigens and attempts to destroy them. \\ When tumor cells grow uncontrollably without control, the general solution of the model becomes $C(t)=\frac{92.6e^{\lambda _{C}t}}{92.6+e^{\lambda _{C}t}}$ neglecting all treatment and activations of immune cells, graphically the solution shown in Fig.(\ref{Fig:3 Graphical prisanation of tumor growth without control}) its linear but an exponential growth. Immune cells can act against tumors in different ways, such as absorbing and presenting tumor antigens, releasing cytokines that activate and recruit other immune cells, or directly killing tumor cells. This act is shown by Fig.(\ref{Fig:3 Graphical presantations of lung tumor growth models}). \begin{center} \begin{figure}[hbt!] \centering %\fbox{% %\tmpframe{ \includegraphics[width = 4.0in]{image1.1.jpg} %} %} \caption{ {\bf Tumor cells population growth with varying intrinsic rate $\lambda _{C}$ denoted by $C(0)$}} \label{Fig:3 Graphical prisanation of tumor growth without control} \end{figure} \end{center} Taking into account the natural death of tumor cells and the activations of immune cells such as natural killer cells along with dendritic cells as presented in Fig.(\ref{Fig:3 Graphical presantations of lung tumor growth models}), such a growth is a straight exponential growth of tumor cells in the lungs reaching a certain limit depending on the reaction of immune cells . Effector T cells are activated and respond to tumor in lungs by migrating to the tumor site, recognizing tumor-associated antigens, and releasing cytotoxins such as granzymes and perforin to induce tumor cell death \cite{vineis2014global}.\\ An effective immune response against tumor infections depends on the activation of cytotoxic T cells that can clear infection by killing tumor infected cells \cite{dwivedi2020perspective}. Proper activation of immune cells depends on professional antigen-presenting cells, such as dendritic cells (DC). Dendritic cell (DC) maturation is inhibited by cigarette smoking, as demonstrated by reduced cell surface expression. Consequently, DCs from animals exposed to cigarette smoke show a reduced capacity to stimulate and activate antigen-specific T-cells in vitro condition; this phenomenon is consistent with reduced antigen-specific T-cell proliferation \cite{saini2020cancer}. %\begin{eqnarray} % \left\{ {\begin{array}{c}\ds C(t)=\frac{d_{D}C_{0}}{\gamma _{C}}\left( {\begin{array}{c}1-e^{\frac{\gamma _{C}t}{2}}\end{array} } \right)~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ % \ds T_{c}(t)=\frac{T_{0}}{4d_{T_{c}}}\left[ {\begin{array}{c}\gamma _{T_{c}}-(\gamma _{T_{c}}-4d_{T_{c}})e^{-d_{T_{c}}t}\end{array} } \right]~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ % \ds D_{T}(t)=\frac{D_{T0}}{2d_{D}}\left[ {\begin{array}{c}\gamma _{D}-(\gamma _{D}-2d_{D})e^{-d_{D}t}\end{array} } \right]~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ % \ds P(t)=\frac{1}{2d_{p}}\left[ {\begin{array}{c}k_{1}\beta _{P}-(k_{1}\beta _{P}-2d_{P}P_{0})e^{-d_{P}t}\end{array} } \right]~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ % \ds F(t)=\frac{d_{P}C_{0}}{\gamma _{C}}\left( {\begin{array}{c}\frac{-\alpha }{r_{F}}-\frac{2\alpha }{2r_{F}+1}\end{array} } \right)e^{\frac{\gamma _{C}t}{2}}-\left( {\begin{array}{c}F_{0}-\frac{\alpha d_{P}C_{0}}{\gamma _{C}r_{F}}+\frac{2\alpha d_{P}C_{0}}{2r_{F}+\gamma _{C}}\end{array} } \right)e^{-r_{F}t}~~~~~~~~~~~\\ % \ds E(t)=\frac{1}{\mu}\left[ {\begin{array}{c}\frac{\gamma _{E}}{4}E_{0}+(1-u)r_{E}F_{0}\end{array} } \right]+e^{-\mu t}\left\{ {\begin{array}{c}E_{0}-\frac{1}{\mu}\left[ {\begin{array}{c}\frac{\gamma _{E}}{4}E_{0}+(1-u)r_{E}F_{0}\end{array} } \right]\end{array} } \right \}\\ % \end{array} } \right.\label{eq49} %\end{eqnarray} \begin{center} \begin{figure}[hbt!] \centering %\fbox{% %\tmpframe{ \includegraphics[width = 3.0in]{image1.2.jpg} \includegraphics[width = 3.0in]{image1.4.jpg} \includegraphics[width = 3.0in]{image1.5.jpg} \includegraphics[width = 3.0in]{image1.6.jpg} %} %} \caption{ {\bf Graphical presentation of tumor cells with control from natural deaths with interaction with Dendritic cells, Antigen peptides and effector T cells}}. \label{Fig:3 Graphical presantations of lung tumor growth models} \end{figure} \end{center} The general solution of the models in Eq. (\ref{eq13}) to Eq. (\ref{eq21}) together with their graphical presentations in Fig. (\ref{Fig:3 Graphical presantations of lung tumor growth models}) shows a stable condition with a predictable state referring to equilibrium state where the cancer is present but not rapidly taking over , this normally happens to a patient with healthy tissues or cells responding quickly to tumor infection in lungs. \\ The quality of the data fitting was also assessed using the root mean square error ($R^{2}$) between the model simulations and between model simulations and actual patient PSA data: \begin{eqnarray} \ds R^{2}&=&1-\frac{SSE}{SST}=1-\frac{\sum _{i}(\widehat{y}_{i}-y_{i})^{2}}{\sum _{i}(\overline{y}_{i}-y_{i})^{2}}\label{eq52} \end{eqnarray} where SSE denotes the sum of squares of errors, SST represents the total sum of squares, $y_{i}$ signifies actual data, $\widehat{y}_{i}$ denotes predicted values, and $\overline{y}_{i}$ stands for the mean of all real data points. The statistic $R^{2}$ determines the degree of fit between the simulation of the model and the actual data, ranging between 0 and 1; higher values indicate better model performance \cite{wang2025optimal}. A negative $R^{2}$ indicates that the predictive model performs worse than a horizontal line representing the mean of the data; essentially, the model's predictions are less accurate than simply guessing the average value for every data point. The individualized parameter values resulting from model fitting are detailed in Table (\ref{T2}) Appendix A shows estimates of growth rate of the tumor cell population of the patients were made and extracted by Zhang (2022) \cite{zhang2022evolution} using time series. As shown in Fig. (\ref{Fig:3 Comparing the carrying capacity}) initially, when the tumor is small and far from the carrying capacity, growth is nearly exponential and the growth rate is the dominant factor \cite{wang2025optimal}. As the tumor approaches carrying capacity ($K_{C}$), the growth rate decelerates and approaches zero, a phenomenon known as saturation. \begin{center} \begin{figure}[hbt!] \centering %\fbox{% %\tmpframe{ \includegraphics[width = 3.0in]{image2.1.jpg} \includegraphics[width = 3.0in]{image2.2.jpg} \includegraphics[width = 3.0in]{image2.3.jpg} \includegraphics[width = 3.0in]{image2.4.jpg} %} %} \caption{ {\bf Comparing the carrying capacity against the tumor growth rate and the tumor death rate for treatment analysis. }} \label{Fig:3 Comparing the carrying capacity} \end{figure} \end{center} According to data in Table (\ref{T3}) of Appendix B, a high resistance rate reduces the activation rate of effector T cells, dendritic cells and estrogen receptors, resulting in a lower administration rate consistent with the resistance rate for a lower concentration. The lower concentration is the reduction of an active pharmaceutical ingredient in plasma, blood, or at the target site over time, usually driven by the metabolic, reversible movement of drugs from the blood to the tissues of the lungs, reduction of plasma concentration or excretory processes \cite{wang2025optimal}. We used the fourth order Runge–Kutta method in numerical simulation to solve the differential equations and to obtain the predicted tumor volume per day shown in Table(\ref{T2}) of Appendix A. This method is well suited for its accuracy and stability in approximating solutions to differential equations, making it a reliable choice for our modeling purpose. \section{Results and Discussions} \subsection{Theoretical analysis} This section theoretically explores the feasibility of the model from Eq.(\ref{eq13}) to Eq.(\ref{eq21}) for a critical mathematical approach to understand tumor immune dynamics, predict treatment outcomes, and optimize therapeutic strategies using stability analysis of the models involving the examination of equilibrium points to determine whether the disease will disappear, persist, or fluctuate over time, directly informs clinical decision making \cite{zhang2022evolution}. For inconvenient mathematical reasons, we then introduced the variables $$(x_{1}, x_{2}, x_{3}, x_{4}, x_{5}, x_{6}, x_{7}, x_{8}, x_{9})=(C, C_{N}, u_{M_{1}}, u_{M_{2}}, E, F, P, D_{T}, T_{c}),$$ \begin{eqnarray} \left \{ {\begin{array}{c}\ds \frac{dx_{1}}{dt}=\gamma _{C}x_{1} \left( {\begin{array}{c} 1-\frac{x_{1}}{K}\end{array} } \right)c_{1}-\beta _{C}\frac{x_{9}^{2}x_{1}}{x_{9}^{2}+M_{T_{c}}} - \beta _D\frac{x_{8}x_{1}}{x_{8}+M_{D}}+\frac{\eta _{C}N_{0}x_{1}}{n}-d_{C}x_{3}x_{1}~~~\\ \ds \frac{dx_{2}}{dt}=\gamma _{C_{N}}x_{2} \left( {\begin{array}{c} 1-\frac{x_{2}}{K}\end{array} } \right)c_{1}-\beta _{C_{N}}\frac{x_{9}^{2}x_{2}}{x_{9}^{2}+M_{T_{c}}} - \beta _D\frac{x_{8}x_{2}}{x_{8}+M_{D}}+\frac{\eta _{C}N_{0}x_{2}}{n}-c_{2}x_{3}x_{2} \\ \ds \frac{dx_{3}}{dt}=\gamma _{u_{M_{1}}}x_{3}\left( {\begin{array}{c} 1-\frac{x_{3}+x_{4}}{K_m}\end{array} } \right)-\alpha _{m_{1}}x_{3}\frac{x_{1}+x_{2}}{x_{1}+x_{2}+M_{u_{M_{1}}}}-d_{m_{1}}x_{3}~~~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{4}}{dt}=\gamma _{u_{M_{2}}}x_{4}\left( {\begin{array}{c} 1-\frac{x_{3}+x_{4}}{K_m}\end{array} } \right)-d_{m_{2}}x_{4}-\alpha _{m_{2}}x_{4}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{5}}{dt}=\lambda _{E} x_{5}(0)\left( {\begin{array}{c}\frac{x_{3}+x_{4}}{x_{3}+x_{4}+M_{E}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{E}}{x_{5}+M_{E}}\end{array} } \right)+(1-u)r_{E}x_{6}-\mu x_{5}~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{6}}{dt} = r_{F}x_{6}\left( {\begin{array}{c}1-\frac{x_{6}}{x_{6}(Max)}\end{array} } \right)-\alpha (x_{1}+x_{2})~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{7}}{dt}=k_{1}\beta _{P}\left( {\begin{array}{c}\frac{x_{9}^{2}}{x_{9}^{2}+M_{T_{c}}}\end{array} } \right)(x_{1}+x_{2})-d_{P} x_{7}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{8}}{dt} = \lambda_{D}x_{8}(0)\left( {\begin{array}{c}\frac{x_{8}+k_{2}(x_{1}+x_{2})}{x_{8}+k_{2}(x_{1}+x_{2})+M_{P}}\end{array} } \right)-d_{D} x_{8}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\\ \ds \frac{dx_{9}}{dt} =\lambda_{T_{c}}x_{9}(0)\left( {\begin{array}{c}\frac{x_{8}}{x_{8}+M_{D}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{T_{c}}}{x_{9}+M_{T_{c}}}\end{array} } \right)-\frac{\beta _{T_{c}}x_{2}x_{9}}{x_{2}+M_{C_{N}}}-d_{T_{c}} x_{9}~~~~~~~~~~~~~~~~~~~~~~~~ \end{array} } \right.\label{eq53} \end{eqnarray} We start by analyzing the existence and stability of the positive steady states of the model. Subsequently, we investigate the feasibility of positive cyclic solutions within the model. This analysis is crucial to assess the potential effectiveness of adaptive cyclic therapy strategies. Also provides a reliable reference point for understanding the highly complex, nonlinear and chaotic dynamics of tumor growth, metabolism, and drug response \cite{uthamacumaran2022review}. Although tumors are generally unstable and non-equilibrium systems, analyzing the theoretical steady state with inputs and outputs balanced to allow mapping of critical thresholds, identifying drug targets, and predict how cells transition to malignant states. \subsubsection{Existence and stability of steady states} We first examine the existence and stability of steady states for Eq.(\ref{eq53}) under constant drug administration dose $\varphi_{0}$ that reduces the immune cell response. letting $\overline{x}^{*}=(x_{1}^{*}, x_{2}^{*}, .......,x_{9}^{*})$ denote the non-negative steady state \cite{wang2025optimal}. To obtain the steady state of the system of Eq.(\ref{eq53}), we need to solve the following equations: \begin{eqnarray} \ds \gamma _{C}x_{1}^{*} \left( {\begin{array}{c} 1-\frac{x_{1}^{*}}{K}\end{array} } \right)c_{1}-\beta _{C}\frac{x_{9}^{*2}x_{1}^{*}}{x_{9}^{*2}+M_{T_{c}}} - \beta _D\frac{x_{8}^{*}x_{1}^{*}}{x_{8}^{*}+M_{D}}+\frac{\eta _{C}N_{0}x_{1}^{*}}{n}-d_{C}x_{3}^{*}x_{1}^{*}=0~~~~~\label{eq54}\\[5pt] \ds \gamma _{C_{N}}x_{2}^{*} \left( {\begin{array}{c} 1-\frac{x_{2}^{*}}{K}\end{array} } \right)c_{1}-\beta _{C_{N}}\frac{x_{9}^{*2}x_{2}^{*}}{x_{9}^{*2}+M_{T_{c}}} - \beta _D\frac{x_{8}^{*}x_{2}^{*}}{x_{8}^{*}+M_{D}}+\frac{\eta _{C}N_{0}x_{2}^{*}}{n}-c_{2}x_{3}^{*}x_{2}^{*} =0~~~~~\label{eq55}\\[5pt] \ds \gamma _{u_{M_{1}}}x_{3}^{*}\left( {\begin{array}{c} 1-\frac{x_{3}^{*}+x_{4}^{*}}{K_m}\end{array} } \right)-\alpha _{m_{1}}x_{3}^{*}\frac{x_{1}^{*}+x_{2}^{*}}{x_{1}^{*}+x_{2}^{*}+M_{u_{M_{1}}}}-d_{m_{1}}x_{3}^{*}=0~~~~~\label{eq56}\\[5pt] \ds \gamma _{u_{M_{1}}}x_{4}^{*}\left( {\begin{array}{c} 1-\frac{x_{3}^{*}+x_{4}^{*}}{K_m}\end{array} } \right)-d_{m_{2}}x_{4}^{*}-\alpha _{m_{2}}x_{4}^{*}=0~~~~~\label{eq57}\\[5pt] \ds \lambda _{E} x_{1}^{*}(0)\left( {\begin{array}{c}\frac{x_{3}^{*}+x_{4}^{*}}{x_{3}^{*}+x_{4}^{*}+M_{E}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{E}}{x_{5}^{*}+M_{E}}\end{array} } \right)+(1-u)r_{E}x_{6}^{*}-\mu x_{5}^{*}=0~~~~~\label{eq58}\\[5pt] \ds r_{F}x_{6}^{*}\left( {\begin{array}{c}1-\frac{x_{6}^{*}}{x_{0}^{*}}\end{array} } \right)-\alpha (x_{1}^{*}+x_{2}^{*})=0~~~~~\label{eq59}\\[5pt] \ds k_{1}\beta _{P}\left( {\begin{array}{c}\frac{x_{9}^{*2}}{x_{9}^{*2}+M_{T_{c}}}\end{array} } \right)(x_{1}^{*}+x_{2}^{*})-d_{P} x_{7}^{*}=0~~~~~\label{eq60}\\[5pt] \ds \lambda_{D}x_{8}^{*}(0)\left( {\begin{array}{c}\frac{x_{8}^{*}+k_{2}(x_{1}^{*}+x_{2}^{*})}{x_{8}^{*}+k_{2}(x_{1}^{*}+x_{2}^{*})+M_{P}}\end{array} } \right)-d_{D} x_{8}^{*}=0~~~~~\label{eq61}\\[5pt] \ds \lambda_{T_{c}}x_{9}^{*}(0)\left( {\begin{array}{c}\frac{x_{8}^{*}}{x_{8}^{*}+M_{D}}\end{array} } \right)\left( {\begin{array}{c}\frac{M_{T_{c}}}{x_{9}^{*}+M_{T_{c}}}\end{array} } \right)-\frac{\beta _{T_{c}}x_{2}^{*}x_{9}^{*}}{x_{2}^{*}+M_{C_{N}}}-d_{T_{c}} x_{9}^{*} =0~~~~~\label{eq62} \end{eqnarray} Starting with $x_{6}^{*}$, $x_{7}^{*}$, $x_{8}^{*}$ and $x_{9}$ with conditions that $x_{1}^{*}=0$, $x_{2}^{*}=0$, then $x_{1}^{*}+x_{2}^{*}=0$ at disease free equilibrium. Then we get $x_{6}^{*}=0$ or $x_{6}^{*}=x_{0}^{*}$ and for $x_{8}^{*}$ as follows \begin{eqnarray} x_{8}^{*}=\frac{\lambda _{D}x_{8}^{*}(0)-M_{P}}{d_{P}} \end{eqnarray} \begin{eqnarray} x_{9}^{*2}=\frac{d_{P}M_{T_{c}}x_{7}^{*}}{k_{1}\beta _{P}(x_{1}^{*}+x_{2}^{*})-d_{P}x_{7}^{*}} \end{eqnarray} From this we obtained $x_{9}^{*}=0$, now its easy to solve \begin{eqnarray} x_{5}^{*}=\frac{-\mu M_{E}+\sqrt{\mu ^{2}M_{E}^{2}-4\mu C}}{2\mu} \label{eq63} \end{eqnarray} Were $C=\frac{\lambda _{E}M_{E}\omega _{M}x_{5}^{*}(0)}{\omega _{M}+M_{E}}$, with \begin{eqnarray} \omega _{M}=x_{3}^{*}+x_{4}^{*}=K_{m} \left[ {\begin{array}{c}1-\frac{d_{m_{2}}+\alpha _{m_{2}}}{\gamma _{u_{M_{1}}}} \end{array} } \right] \label{eq64} \end{eqnarray} Therefore, from the given conditions of $x_{1}$ and $x_{2}$ we get positive values as in Eq.(\ref{eq64}) that is $x_{3}^{*}+x_{4}^{*}\geq 1$ and Eq.(\ref{eq63}) $\sqrt{\mu ^{2}M_{E}^{2}-4\mu C}\geq \mu M_{E}$, there exists a unique positive solution, hence we represent the steady state $\overline{x}^{*}=(x_{1}^{*}, x_{2}^{*}, .......,x_{9}^{*})$ of the system in Eq.(\ref{eq53}) using cancer cell numbers $E=(x_{1}^{*}, x_{2}^{*}, x_{3}^{*})$ \begin{theorem} Consider the system in Eq.(\ref{eq53}) with constant drug dose rate $\varphi_{0}>0$ increasing the response rate of immune cells, leading to tumor cells $x_{1}^{*}$ and $x_{2}^{*}$ becoming greater than zero. \begin{itemize} \item [(1)] There always exists a zero steady state $E_{0}=(0,0,0)$ \item [(2)] If $\frac{\alpha _{m_{3}} }{d_{m_{3}}}<1$ and also that $\frac{d_{m_{2}}}{\gamma _{u_{M_{2}}}}<\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}}\epsilon \Bbb{R}^{+}$, there exists a steady state $$E_{1}= \left\{ {\begin{array}{c}0~~,0~~, \frac{K_{m}}{a-1}\left[ {\begin{array}{c}\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}}-\frac{d_{m_{2}}}{\gamma _{u_{M_{2}}}}\left( {\begin{array}{c}1-\frac{\alpha _{m_{3}}}{d_{m_{3}}}\end{array} } \right)\end{array} } \right]~~~\end{array} } \right\}$$ \item [(3)] Define the functions $f_{2}:\Bbb{R}\mapsto\Bbb{R}$ and $f_{1}:\Bbb{R}\mapsto\Bbb{R}$ as \begin{eqnarray} f_{2}(x_{3})=K+c_{6}(c_{7}-c_{2}x_{3}) \end{eqnarray} and \begin{eqnarray} f_{1}(x_{3})=K+c_{3}(c_{4}-c_{5}x_{3}) \end{eqnarray} were \begin{eqnarray} \ds &&c_{1}=1+r_{m_{2}}x_{4},~~c_{3}=\frac{K}{\gamma _{C}c_{1}},~~c_{4}=\frac{\eta _{C}N_{0}}{n}\cr && c_{5}=d_{C},~~c_{6}=\frac{K}{\gamma _{C_{N}}c_{1}},~~c_{7}=\frac{\eta _{C_{N}}N_{0}}{n},~~ c_{2}=e_{C_{N}}+d_{C_{N}} \end{eqnarray} The positive steady state $E_{2}=(x_{1}^{*}, x_{2}^{*}, x_{3}^{*})$ only exists if and only if there exists $x_{1}>1$ such that $$f_{2}(x_{3}^{*})>0$$ and satisfies the equation $$x_{1}=f_{2}\circ \widehat{f_{3}}(x_{1}^{*})$$ were $\circ$ represent the composition of the function, and $\widehat{f}:\Bbb{R}\mapsto \Bbb{R}^{2}$ is defined as $$\widehat{f}(x_{1})=(x_{1}, x_{2})$$ Here, $f_{3}(x_{1})$ is defined by Eq(\ref{eq54}). In this case, the steady state is given by $$E_{2}=(x_{1}^{*}, f(x_{1}^{*}, \widehat{f_{3}}(x_{1})))$$ \end{itemize} \label{theoem1.2} \end{theorem} \begin{proof} (1) From Eq.(\ref{eq54}) to Eq.(\ref{eq55}), it is easy to have an equilibrium state $E_{0}=(0, 0, 0)$. (2) From Eq(\ref{eq55}), when $x_{1}^{*}=x_{2}^{*}=0$ neglects the terms with $x_{3}^{*}, x_{8}^{*}, x_{9}^{*}$ with $c_{1}=1$, we have \begin{eqnarray} \ds && \gamma _{u_{M_{1}}}x_{3} \left( {\begin{array}{c}1-\frac{x_{3}+x_{4}}{K_{m}}\end{array} } \right)-d_{m_{3}}x_{3}=0 \cr && \ds \gamma _{u_{M_{2}}}x_{4}\left( {\begin{array}{c}1-\frac{x_{3}+x_{4}}{K_{m}}\end{array} } \right)-d_{m_{3}}x_{4}-\alpha _{m_{3}}x_{4}=0\label{eq65} \end{eqnarray} Then we obtain the two equations of $x_{3}+x_{4}$ from Eq.(\ref{eq65}) with the coefficients of $x_{3}$ as $a$ so that $a>1$ for a steady state condition, \begin{eqnarray} \ds && ax_{3}+x_{4}=K_{m} \left[ {\begin{array}{c}1-\frac{d_{m_{3}}}{\gamma _{u_{M_{2}}}} \left( {\begin{array}{c}1-\frac{\alpha _{m_{3}}}{d_{m_{3}}}\end{array} } \right)\end{array} } \right]\cr && x_{3}+x_{4}=K_{m}\left[ {\begin{array}{c}1-\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}}\end{array} } \right]\label{eq1} \end{eqnarray} Solving the equations in (\ref{eq1}) simultaneously gives the steady state \begin{eqnarray} E_{1}= \left\{ {\begin{array}{c}0~~,0~~, \frac{K_{m}}{a-1}\left[ {\begin{array}{c}\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}}-\frac{d_{m_{2}}}{\gamma _{u_{M_{2}}}}\left( {\begin{array}{c}1-\frac{\alpha _{m_{3}}}{d_{m_{3}}}\end{array} } \right)\end{array} } \right]~~~\end{array} } \right\}\label{eq66} \end{eqnarray} (3) Considering a steady state $(x_{1}^{*}, x_{2}^{*})$ with $x_{1}^{*}>0$, then $f_{1}(x_{3})f^{-1}_{2}(x_{3})=1$ and $f_{3}(x_{1})f_{3}^{-1}(x_{2})=1$ were $x_{1}=x_{2}$ thus the positive steady state $E_{2}=(x_{1}^{*},x_{2}^{*},x_{3}^{*})$ exists if and only if the equation $x_{1}^{*}=f_{2}\circ \widehat{f}_{3}(x_{1}^{*})$ has a positive solution. From the proof of Theorem (\ref{theoem1.2}), we reduce the existence of positive $x_{1}$ and $x_{2}$ steady states into a problem of finding the positive solution of a single equation. Many numerical methods can be applied to solve a single equation. Nevertheless, obtaining the explicit condition for the model parameters is not easy \cite{wang2025optimal}. Clinically, steady state $E_{1}$ means the dominant situation of resistant cells, while steady state $E_{2}$ means the coexistence of sensitive and resistant cancer cells \end{proof} \begin{theorem} Consider the system in Eq.(\ref{eq53}) with a high immune response as a result of drug administration $D_{j}(t)$ with rate $\varphi \geq 1$, the zero steady state $E_{0}$ is unstable. Furthermore,the steady-state $E_{1}$ is locally asymptotically stable if and only if \begin{eqnarray} \gamma _{C}(1+r_{m_{2}}f_{4}(0))-\frac{\beta _{C}f_{9}(0)^{2}}{f_{9}(0)^{2}+M_{T_{c}}}-\frac{\beta _{D}f_{8}(0)}{f_{8}(0)+M_{D}}+\frac{\eta _{C}N_{0}}{n}-d_{C}f_{3}(0)<0\label{eq67} \end{eqnarray} and \begin{eqnarray} \gamma _{C_{N}} (1+r_{m_{2}}f_{4}(0))-\frac{\beta _{C_{N}}f_{9}(0)^{2}}{f_{9}(0)^{2}+M_{T_{c}}}-\frac{\beta _{D}f_{8}(0)}{f_{8}(0)+M_{D}}+\frac{\eta _{C_{N}}N_{0}}{n}-c_{2}f_{3}(0)<0\label{eq68} \end{eqnarray} \end{theorem} \begin{proof} Firstly, we state that the linear stability of steady states is controlled by the eigenvalues of the Jacobian matrix associated with the system in Eq.(\ref{eq53}) $$J(x_{1}^{*}, x_{2}^{*},x_{3}^{*},x_{4}^{*},x_{5}^{*},x_{6}^{*},x_{7}^{*},x_{8}^{*},x_{9}^{*})= \left. {\begin{array}{c c c c c c c c c} \left( {\begin{array}{c c c c c c c c c} J_{11} & 0 & J_{13} & 0 & 0 & 0 & 0 & J_{18} & J_{19}\\ 0 & J_{22} & J_{23} & 0 & 0 & 0 & 0 & J_{28} & J_{29}\\ J_{31} & J_{32} & J_{33} & J_{34} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & J_{43} & J_{44} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & J_{53} & J_{54} & J_{55} & J_{56} & 0 & 0 & 0 \\ J_{61} & J_{62} & 0 & 0 & 0 & J_{66} & 0 & 0 & 0 \\ J_{71} & J_{72} & 0 & 0 & 0 & 0 & J_{77} & 0 & J_{79}\\ J_{81} & J_{82} & 0 & 0 & 0 & 0 & 0 & J_{88} & 0 \\ 0 & J_{92} & 0 & 0 & 0 & 0 & 0 & J_{98} & J_{99}\\ \end{array} } \right)\end{array} } \right|_{x^{*}=\overline{x}^{*}}$$ where $$J_{11}=\gamma _{C} \left( {\begin{array}{c}1-\frac{2x_{1}}{K}\end{array} } \right)c_{1}-\frac{\beta _{C}x_{9^{2}}}{x_{9}^{2}+M_{T_{c}}}-\frac{\beta _{D}x_{8}}{x_{8}+M_{D}}+\frac{\eta _{C}N_{0}}{n}-d_{C}x_{3}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{13}=-d_{C}x_{1}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{18}=-\beta _{D}\frac{x_{1}}{x_{8}+M_{D}}-\beta _{D}\frac{x_{1}x_{8}}{(x_{8}+M_{D})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{19}=-2\beta _{C}\frac{x_{9}x_{1}}{x_{9}^{2}+M_{T_{c}}}-2\beta _{C}\frac{x_{9}^{2}x_{1}}{(x_{9}^{2}+M_{T_{C}})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{22}=\gamma _{C_{N}} \left( {\begin{array}{c}1-\frac{2x_{2}}{K}\end{array} } \right)c_{1}-\frac{\beta _{C}x_{9}^{2}}{x_{9}^{2}+M_{T_{c}}}-\frac{\beta _{D}x_{8}}{x_{8}+M_{D}}+\frac{\eta _{C}N_{0}}{n}-c_{2}x_{3}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{23}= -c_{2}x_{2}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{28}=-\beta _{D}\frac{x_{2}}{x_{8}+M_{D}}-\beta _{D}\frac{x_{2}x_{8}}{(x_{8}+M_{D})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{29}=-2\beta _{C_{N}}\frac{x_{9}x_{2}}{x_{9}^{2}+M_{T_{c}}}-2\beta _{C_{N}}\frac{x_{9}^{2}x_{2}}{(x_{9}^{2}+M_{T_{C}})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{31}=\frac{-1}{\alpha _{m_{1}}(x_{1}+x_{2})}-\frac{x_{1}+x_{2}+M_{u_{M_{1}}}}{\alpha _{m_{1}}x_{3}(x_{1}+x_{2})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{32}=\frac{-1}{\alpha _{m_{1}}(x_{1}+x_{2})}-\frac{x_{1}+x_{2}+M_{u_{M_{1}}}}{\alpha _{m_{1}}x_{3}(x_{1}+x_{2})^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{33}=\gamma _{u_{M_{1}}}\left( {\begin{array}{c}1-\frac{2x_{3}}{K_{m}}\end{array} } \right)-d_{m_{1}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{34}=\frac{-\gamma _{u_{M_{1}}}x_{3}}{K_{m}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{43}=\frac{-\gamma _{u_{M_{2}}}x_{4}}{K_{m}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{44}=\gamma _{u_{M_{2}}} \left( {\begin{array}{c}1-\frac{x_{3}+2x_{4}}{K_{m}}\end{array} } \right)-d_{m_{2}}-\alpha _{m_{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{53}=\frac{x_{5}+M_{E}}{\gamma _{E}x_{5}(0)M_{E}}\left[ {\begin{array}{c}1+\frac{M_{E}+1}{(x_{3}+x_{4})^{2}}\end{array} } \right]~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{54}=\frac{\gamma _{E}x_{5}(0)M_{E}}{x_{5}+M_{E}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{55}=\frac{1}{\gamma _{E}x_{5}(0)M_{E}}\left[ {\begin{array}{c}1+\frac{1}{x_{3}+x_{4}}\end{array} } \right]~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{56}=(1-u)r_{E}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{61}=-\alpha~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{62}=-\alpha~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{66}=r_{F}\left( {\begin{array}{c}1-\frac{2x_{6}}{x_{0}}\end{array} } \right)~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{71}=k_{1}\beta _{P}\left( {\begin{array}{c}\frac{x_{9}^{2}}{x_{9}^{2}+M_{T_{c}}}\end{array} } \right)~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{72}=k_{1}\beta _{P}\left( {\begin{array}{c}\frac{x_{9}^{2}}{x_{9}^{2}+M_{T_{c}}}\end{array} } \right)~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{77}=d_{P}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{79}=\frac{2(2x_{9}^{2}+M_{T_{c}})}{k_{1}\beta _{P}x_{9}^{3}(x_{1}+x_{2})}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{81}=\frac{\lambda _{D}k_{2}x_{8}(0)}{x_{8}+k_{2}(x_{1}+x_{2})+M_{P}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{82}=\frac{\lambda _{D}k_{2}x_{8}(0)}{x_{8}+k_{2}(x_{1}+x_{2})+M_{P}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{88}=\frac{\lambda _{D}x_{8}(0)}{x_{8}+k_{2}(x_{1}+x_{2})+M_{P}}-d_{D}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{92}=\frac{2x_{2}+M_{C_{N}}}{\beta _{T_{C}}x_{9}x_{2}^{2}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{98}=\frac{\lambda _{T_{c}}x_{9}(0)}{x_{8}+M_{D}}\left( {\begin{array}{c}\frac{M_{{T}_{c}}}{x_{9}+M_{T_{c}}}\end{array} } \right)~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ $$J_{99}=-\lambda _{T_{c}}x_{9}(0)\left( {\begin{array}{c}\frac{x_{8}}{x_{8}+M_{D}}\end{array} } \right)\frac{M_{T_{c}}}{(x_{9}+M_{T_{c}})^{2}}-\frac{\beta _{T_{c}}x_{2}}{x_{2}+M_{C_{N}}}-d_{T_{c}}~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~$$ For a steady state with $x_{1}^{*}=0$, $x_{2}^{*}=0$ assuming that $\alpha =0$ is a disease free equilibrium state where a disease is completely absent, that is all cells infected, exposed, and recovered compartments are zero, leaving only the susceptible population, then $$ \Bbb{J}= \left( {\begin{array}{c c c c c c c c c} J_{11} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & J_{22} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & J_{33} & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & J_{44} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 &0 & 0 & J_{55} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & J_{66} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & J_{77} & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & J_{88} & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & J_{99}\\ \end{array} } \right) $$ Then we have $Det(\Bbb{J}-\lambda \Bbb{I})=0$ were $I$ the identity matrix to get for eigen values $$(J_{11}-\lambda _{1})(J_{22}-\lambda _{2})(J_{33}-\lambda _{3})(J_{44}-\lambda _{4})(J_{55}-\lambda _{5})(J_{66}-\lambda _{6})(J_{77}-\lambda _{7})(J_{88}-\lambda _{8})(J_{99}-\lambda _{9})=0$$ Hence we get the eigen value $ \overline{\lambda _{i}}$for $\overline{\lambda _{i}}=(\lambda _{1}, ....., \lambda _{9})$ and $\overline{J_{ij}}=(J_{11}, .......,J_{99})$ with $(J_{22}, .....j_{99})<0$ then $\overline{\lambda _{i}}=\overline{J_{ij}}$ in the steady state $E_{0} = (0, 0, 0)$ therefore $E_{0}$ is unstable. When $\frac{\eta _{C}N_{0}}{\gamma_{C_{N}n}}\leq 0$ and consider the steady state $\frac{\alpha _{m_{3}} }{d_{m_{3}}}<1$ and also that $\frac{d_{m_{2}}}{\gamma _{u_{M_{2}}}}<\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}} \epsilon \Bbb{R}^{+}$ considering the steady state $E_{1}=(0,0,x_{3}^{*})$ we have $$x_{3}^{*}=\frac{K_{m}}{a-1}\left[ {\begin{array}{c}\frac{d_{m_{1}}}{\gamma _{u_{M_{1}}}}-\frac{d_{m_{2}}}{\gamma _{u_{M_{2}}}}\left( {\begin{array}{c}1-\frac{\alpha _{m_{3}}}{d_{m_{3}}}\end{array} } \right)\end{array} } \right]$$ and hence $$\lambda _{3}=\gamma _{u_{M_{1}}}\left( {\begin{array}{c}1-\frac{2x_{3}^{*}}{K_{m}}\end{array} } \right)-d_{m_{1}}<0~~~~~~~~~~$$ In steady state $E_{1}$, we have $x_{3}^{*}=f_{3}(0),x_{4}^{*}=f_{4}(0), x_{8}^{*}=f_{8}(0)$ and $x_{9}^{*}=f_{9}(0)$ with \begin{eqnarray} \ds && \lambda _{1}= \gamma _{C}(1+r_{m_{2}}f_{4}(0))-\frac{\beta _{C}f_{9}(0)^{2}}{f_{9}(0)^{2}+M_{T_{c}}}-\frac{\beta _{D}f_{8}(0)}{f_{8}(0)+M_{D}}+\frac{\eta _{C}N_{0}}{n}-d_{C}f_{3}(0)<0\cr && \ds \lambda _{2}= \gamma _{C_{N}} (1+r_{m_{2}}f_{4}(0))-\frac{\beta _{C_{N}}f_{9}(0)^{2}}{f_{9}(0)^{2}+M_{T_{c}}}-\frac{\beta _{D}f_{8}(0)}{f_{8}(0)+M_{D}}+\frac{\eta _{C_{N}}N_{0}}{n}-c_{2}f_{3}(0)<0\label{eq69} \end{eqnarray} Thus, the steady-state $E_{1}$ is locally asymptotically stable if and only if $\lambda _{1} < 0$. \end{proof} \subsection{Optimizing Chemotherapy Medications} Lung cancer chemotherapy contains many side effects that are also part of the treatment process, resulting in the reduction of effector T cell number that increases the risk of cancer effects. This involves having a low number of cells that help fight infection, not only effector T cells but other cells such as dendritic cells. According to Xiong at el \cite{xiong2025mathematical} when therapies have serious side effects, the objective function is composed of two parts denoted by $\phi _{1}(S)$ which is the tumor risk function, the tumor clone risk function $\phi _{2}(S)$ and $\phi _{3}(S)$ the treatment risk function that is $\phi (S)=\phi _{1}(S)+\phi _{2}(S)+\phi _{3}(S)$ then, \begin{eqnarray} \phi _{1}(S)=\omega _{1}\frac{C_{tfinal}-C_{min}}{C_{max}-C_{min}}\label{eq70} \end{eqnarray} \begin{eqnarray} \phi _{2}(S)=\omega _{2}\frac{C_{(N; tfinal)}-C_{(N; min)}}{C_{(N; max)}-C_{(N; min)}}\label{eq71} \end{eqnarray} where $\omega _{1}$ is the weight of the tumor risk, $\omega _{2}$ is the weight of the risk of tumor clones, $C_{max}$ is the final size of the untreated tumor, $C_{(N; max)}$ is the size of the final size of the untreated tumor clones, $C_{min}$ is the minimum final tumor volume subject to drug constraints and $C_{(N; min)}$ is the minimum final tumor clone volume subject to drug constraints. Then for $\phi _{3}$ we generate it from dose toxicity curve illustrating that toxicity increases with higher concentrations of a drug. The relationship is often mathematically described using a logistic model then \begin{eqnarray} \phi _{3}(S)=\omega _{3}r_{s}\label{72} \end{eqnarray} were $\omega _{3}$ is the weight of tumor and $r_{s}$ is the sigmoidal dose response equation given by \begin{eqnarray} r_{s}=r_{a}+\frac{r_{b}-r_{a}}{1+10^{\left( {\begin{array}{c}\frac{log(r_{1})}{r_{0}}-r_{t}\end{array} } \right) h}}\label{73} \end{eqnarray} with $r_{b}$ the maximum toxicity plateaus of the curve, $r_{a}$ minimum toxicity plateaus of the curve, $r_{1}=EC_{50}$ doses that produce $50\%$ of the maximum toxic effects, $r_{0}$ is the dose, $r_{t}$ is the log of the dose or concentration ($\log _{10}M$) and $h$ is the hill slope which defines the steepness of the curve \cite{manukonda2025dose, greim2018introduction}. The dose equation is used to calculate efficacy and potency in in vitro cell assays. It shows that increasing drug concentration leads to a steeper decrease in cell viability until a maximum effect is reached. According to to Xiong at el \cite{xiong2025mathematical} from the S shaped dose toxicity curve in Fig.(\ref{Fig:4 dose toxicity curve}), once the plateau is reached, increasing the dose of drugs will not enhance efficacy. \begin{center} \begin{figure}[hbt!] \centering \fbox{% %\tmpframe{ \includegraphics[width = 4.0in]{image2.5.jpg} %} } \caption{ {\bf The S shaped dose toxicity curve \cite{manukonda2025dose} }} \label{Fig:4 dose toxicity curve} \end{figure} \end{center} To characterize the phenomenon, we introduce an indicator \begin{eqnarray} r=\frac{r_{max}.r_{t}^{n}}{r_{1}^{n}+r_{t}^{n}}\label{eq72} \end{eqnarray} with $r_{max}$ as the maximal response the drug can produce (efficacy) and $n$ as the slope coefficient of the hill that describes the steepness of the curve. When a Hill coefficient $n>1$ suggests positive cooperative in binding which increases affinity for further binding, while $n=1$ suggests non-cooperative, simple Michaelis-Menten kinetics. The dose-response curve act as a typically semi-logarithmic with the X-axis being the logarithm of the concentration or dose and the Y axis being response, usually measured as a percentage of some baseline value \cite{tallarida2012dose}. \begin{remark} When optimizing a drug, since it is considered to be safe, the optimization function $\phi (S)=\phi _{1}(S)+\phi _{2}(S)$ because both $C_{max}$ and $C_{min}$ are constants in Eq.(\ref{eq70}) and also for $c_{(N; max)}$ and $C_{(N; min)}$ both are constants in Eq.(\ref{eq71}) then the objective function $\phi (S)=C(t_{final})+C_{N}(t_{final})$. However, when optimizing the dose of the drug, the objectives function inevitably results in an unbounded increase in the drug dose. Hence the presents of the drug efficacy plateau allows Eq(\ref{eq72}) to be used to constrain the unlimited increase in the drug dose. \end{remark} %\subsection{Implications of Single Drug and Immune Response} %\subsection{Effects Combination Therapy} %\subsection{Effects of Delaying Second Drug} %\subsection{Treatment interruptions} \section{Conclusion and Discussion} The continuous development of computer technology and mathematical models with the application of mathematical models in medicine has become more and more widespread. Biomedical models based on differential equations, medical models based on statistics, medical models based on machine learning, and medical models based on network science have become indispensable tools in the medical field \cite{Sun}. In this paper, a review and outlook of recent medical problems based on mathematical models are presented from the perspective of mathematical models. Biomedical models are widely used for the early diagnosis of diseases, drug delivery \cite{Idrees2021}. Most lung carcinomas are diagnosed at an advanced stage, conferring a poor prognosis. The need to diagnose lung cancer at an early and potentially curable stage is therefore obvious, and this has inspired the adoption of lung cancer screening in patients at high risk for lung cancer \cite{spiro2010lung}.\\ Patients who develop lung cancer have been smokers and have smoking-related damage to the heart and lungs, making aggressive surgical or multi-modal therapies less viable options. Since drug therapy represents a promising and safe approach in cancer treatment with little side effects compared to other treatment options, the approach diverges from traditional continuous therapy strategies by dynamically adjusting treatment regimens based on patient responses in real time \cite{xiong2025mathematical}. In this work, a mathematical model is constructed to describe the interaction between tumor cells and immune cells in the lungs under the combination chemotherapy drugs that activate the immune cells to fight against infection. The mathematical model was used to predict the volume of tumor cells per day using the Runge Kutta method and was also used to optimize the chemotherapy medication model, reducing tumor size and drug related side effects \cite{manukonda2025dose}. \\ Systematic modeling and quantitative analysis are essential for investigating the dynamic changes within the tumor in lungs. By developing detailed tumor-immune regulatory network models, researchers can quantitatively represent interactions between tumors and the immune system, which aids in identifying potential immune biomarkers predictive of tumor behavior \cite{xiong2025mathematical}. Quantitative metrics derived from these models provide theoretical foundations for understanding cancer immune editing and classifying cancer immune phenotypes. Such metrics not only shed light on tumor-immune system evolution but also facilitate cancer sub-typing \cite{ding2025characterising}. Ultimately, systematic modeling and quantitative analysis offer novel perspectives for cancer therapy, significantly supporting individualized treatment plans for patients \cite{manukonda2025dose}. Cancer is inherently heterogeneous, exhibiting variations in genetics, cellular phenotypes, and interactions with immune and stromal components. Traditional experimental methods alone often fail to capture these multi-scale complexities \cite{azizi2024mathematical}.\\ Our model reduce the complexity of tumor immune interaction with the procedure that the independent of both the shape of the phenotypic distribution and the functional form of the terms that characterize adaptation. This generality and flexibility lends our analysis suitable for adaptation in a wide range of contexts, both within mathematical oncology and more broadly \cite{villa2025reducing}. We expect many of the advantages conferred by the reduced model to become even more pertinent in high-dimensional phenotype spaces, particularly in the context of mechanistic interpretation of correspondingly high-dimensional multi-omics data. The necessity to impose a system closure yields an approximation of the underlying dynamics. However, we highlight that the presented approach can be applied up to an arbitrary order \cite{ding2025characterising}.\\ Optimizing the choice of drugs for different stages of the disease is important when developing a medication plan. Additionally, the tumor cells may become resistant to drugs, which makes it difficult to control the growth of tumor cells \cite{xiong2025mathematical}. The resistance may be induced by the genetic variation in tumor cells, the influence of tumor microenvironment, and other factors. Investigating therapeutic regimens for drug-resistant tumors is a focus of our future research \cite{zhang2022evolution}. Cancer treatment resistance arises from a confluence of intrinsic genetic alterations, tumor microenvironment influences, and other cellular and systemic factors. \subsubsection*{\bf Data Availability Statement} No data was used for the research described in the article. \subsubsection*{\bf Conflicts of Interest} The authors declare that they have no conflicts of interest. \subsubsection*{\bf Author`s Contribution} Author contributed in formulating lung cancer models helping to predict the spread of infection to normal cells by simplifying complex system translating real world phenomena into equations allowing experts to gain deeper insights into their functioning. % --- Appendix starts here --- \newpage \section*{APPENDIX A} \begin{table}[!h] \caption{Estimates of sensitive, tumor growth rates from patient data in both the contemporaneous and adaptive therapy cohorts.} \label{T2} \small\small \begin{center} \begin{tabular}{l l l l l l l l l} \hline \\ Patient & Extracted & Carrying & Clarence rate & Natural & Saturation & Saturation & Tumor & Predicted \\ identifier & growth rate & capacity & Maximum & death rate & coefficient & coefficient & Volume & Tumor\\ & ($\gamma _{C}$) & ($K_{C}$) & ($\beta _{C}$) & ($d_{C}$) & ($M_{D}$) & ($M_{T_{c}}$)& per day & Volume\\ & day$^{-1}$ & $mm^{3}$ & day$^{-1}$ & day$^{-1}$ & & & ($\frac{dC}{dt}$) & per day \\ \hline\\ $C001$ & $0.0022$ & $0.002626$ & $8.1\times 10^{-2}$ & $0.016626$ & $288.702$ &$2240$ & $38.848$ & $47.553$\\ $C002$ & $0.0031$ & $0.003230$ & $8.1\times 10^{-2}$ & $0.033940$ & $141.423$ & $2240$& $38.713$ & $41.500$\\ $C003$ & $0.0100$ & $0.010547$ & $8.1\times 10^{-2}$ & $0.440302$ & $10.9016$ & $2240$& $38.329$ & $40.052$\\ $C004$ & $0.0196$ & $0.024121$ & $8.1\times 10^{-2}$ & $2.352818$ & $2.04010$ & $2240$& $36.579$ & $46.248$\\ $C005$ & $0.0173$ & $0.020588$ & $8.1\times 10^{-2}$ & $1.767238$ & $2.71610$ & $2240$& $37.115$ & $44.627$\\ $C006$ & $0.0130$ & $0.014115$ & $8.1\times 10^{-2}$ & $0.822589$ & $5.83520$ & $2240$& $37.980$ & $41.407$\\ $C007$ & $0.0062$ & $0.006384$ & $8.1\times 10^{-2}$ & $0.149634$ & $32.0782$ & $2240$& $38.590$ & $40.900$\\ $C008$ & $0.0047$ & $0.004808$ & $8.1\times 10^{-2}$ & $0.082004$ & $58.5335$ & $2240$& $38.649$ & $39.964$\\ $C009$ & $0.0107$ & $0.011351$ & $8.1\times 10^{-2}$ & $0.515913$ & $9.30390$ & $2240$& $38.260$ & $41.762$\\ $C010$ & $0.0071$ & $0.007345$ & $8.1\times 10^{-2}$ & $0.201961$ & $23.7669$ & $2240$& $38.544$ & $41.039$\\ $C011$ & $0.0123$ & $0.013251$ & $8.1\times 10^{-2}$ & $0.719190$ & $6.67418$ & $2240$& $38.075$ & $42.198$\\ $C012$ & $0.0075$ & $0.007777$ & $8.1\times 10^{-2}$ & $0.228278$ & $21.0270$ & $2240$& $38.520$ & $41.109$\\ $C013$ & $0.0109$ & $0.011583$ & $8.1\times 10^{-2}$ & $0.538948$ & $8.90623$ & $2240$& $38.239$ & $41.810$\\ $C014$ & $0.0142$ & $0.015654$ & $8.1\times 10^{-2}$ & $1.022376$ & $4.69495$ & $2240$& $37.797$ & $31.674$\\ $P1001$ & $0.0317$ & $0.193273$ & $8.1\times 10^{-2}$ & $9.675990$ & $0.49607$& $2240$& $29.981$ & $29.981$\\ $P1002$ & $0.0245$ & $0.037463$ & $8.1\times 10^{-2}$ & $4.394687$ & $1.09222$& $2240$& $34.719$ & $54.049$\\ $P1003$ & $0.0105$ & $0.011119$ & $8.1\times 10^{-2}$ & $0.493523$ & $9.72599$ & $2240$& $38.280$ & $31.876$\\ $P1004$ & $0.0071$ & $0.007345$ & $8.1\times 10^{-2}$ & $0.201961$ & $23.7670$ & $2240$& $38.544$ & $40.851$\\ $P1006$ & $0.0216$ & $0.028466$ & $8.1\times 10^{-2}$ & $3.070980$ & $1.56302$ & $2240$& $35.923$ & $48.304$\\ $P1007$ & $0.0118$ & $0.012646$ & $8.1\times 10^{-2}$ & $0.650886$ & $7.37457$ & $2240$& $38.137$ & $41.846$\\ $P1011$ & $0.0304$ & $0.106172$ & $8.1\times 10^{-2}$ & $8.460746$ & $0.56733$ & $2240$& $31.062$ & $31.062$\\ $P1012$ & $0.0160$ & $0.018138$ & $8.1\times 10^{-2}$ & $1.381082$ & $3.47554$ & $2240$& $37.469$ & $43.446$\\ $P1014$ & $0.0109$ & $0.011583$ & $8.1\times 10^{-2}$ & $0.538948$ & $8.90623$ & $2240$& $38.239$ & $41.610$\\ $P1015$ & $0.0419$ & $-0.02242$ & $8.1\times 10^{-2}$ & $25.47936$ & $0.18839$ & $2240$& $16.346$ & $16.346$\\ $P1016$ & $0.0106$ & $0.011235$ & $8.1\times 10^{-2}$ & $0.504638$ & $9.51177$ & $2240$& $38.270$ & $41.588$\\ $P1017$ & $0.0068$ & $0.007024$ & $8.1\times 10^{-2}$ & $0.183480$ & $26.1608$ & $2240$& $38.560$ & $40.807$\\ $P1018$ & $0.0189$ & $0.022811$ & $8.1\times 10^{-2}$ & $2.133840$ & $2.24947$ & $2240$& $36.779$& $45.357$\\ \hline \end{tabular} \end{center} \end{table} %{\bf NB:}\\ %\fbox{ %$K_{C}=\frac{\gamma _{C}C[1+e^{(\gamma _{C}-d_{C})}]}{C-C_{0}e^{(\gamma _{C}-d_{C})t}}$ %}\\\\ %\fbox{$d_{C}=\frac{C(t)\gamma _{C}}{C_{0}(1-e^{\frac{\gamma _{C}t}{2}})}$}\\\\ Conditions: $t=120$ days, $C_{0}=0.9893~~mm^{3}$ and then $C(120)=53~~mm^{3}$, maximum clearance rate by effector T cells, $\beta _{C}=8.1\times 10^{-2}$ \newpage \section*{APPENDIX B} \begin{table}[!h] \caption{Estimates of drug resistance with different blood densities from normal to tumor infected blood} \label{T3} \small\small \begin{center} \begin{tabular}{l l l l l l l } \hline \\ Dosage& Mass& Mass & Resistance rate& Resistance rate& Resistance rate& Resistance rate\\ $j$ & $mg$ & $kg$ & $\rho =1060kgm^{-3}$&$\rho =2120kgm^{-3}$& $\rho =3180kgm^{-3}$& $\rho =4240kgm^{-3}$\\ \hline $1$ & $500$ & $0.0005$& $0.004336$ & $0.008673$ & $0.013009$ & $0.017345$\\ $2$ & $450$ & $0.00045$& $0.003903$ & $0.007805$ & $0.011708$ & $0.015611$\\ $3$ & $400$ & $0.00040$& $0.003469$ & $0.006938$ & $0.010407$ & $0.013876$\\ $4$ & $350$ & $0.00035$& $0.003035$ & $0.006071$ & $0.009106$ & $0.012142$\\ $5$ & $300$ & $0.00030$& $0.002602$ & $0.005204$ & $0.007805$ & $0.010407$\\ $6$ & $250$ & $0.00025$& $0.002168$ & $0.004336$ & $0.006505$ & $0.008673$\\ $7$ & $200$ & $0.00020$& $0.001735$ & $0.003469$ & $0.005204$ & $0.006938$\\ $8$ & $150$ & $0.00015$& $0.001301$ & $0.002602$ & $0.003903$ & $0.005204$\\ $9$ & $100$ & $0.00010$& $0.000867$ & $0.001735$ & $0.002602$ & $0.003469$\\ $10$ & $50$ & $0.00005$& $0.000434$ & $0.000867$ & $0.001301$ & $0.001735$\\ \hline \end{tabular} \end{center} \end{table} \fbox{$a_{r}=\frac{\rho v LM}{\mu \Delta t}f(\Pi _{s})$} \begin{table}[!h] \caption{Estimates of concentration of the administered drugs with different masses for analysis} \label{T4} \small\small \begin{center} \begin{tabular}{l l l l l l l } \hline \\ Time($t$) & $D_{j}(t)$ & $D_{j}(t)$ & $D_{j}(t)$ & $D_{j}(t)$ & $D_{j}(t)$ & $D_{j}(t)$\\ ($s$) & $500mg$ & $450mg$ & $400mg$ & $350mg$ & $300mg$ & $250mg$\\ \hline $0$ & $500$ & $450$ & $400$ & $350$ & $300$ & $250$ \\ $1$ & $495.6823863$ & $446.1141477$ & $396.545909$ & $346.9776704$ & $297.4094318$ & $247.8411932$\\ $2$ & $491.4020562$ & $442.2618506$ & $393.1216449$ & $343.9814393$ & $294.8412337$ & $245.7010281$\\ $3$ & $487.1586877$ & $438.4428189$ & $389.7269502$ & $341.0110814$ & $292.2952126$ & $243.5793438$\\ $4$ & $482.9519616$ & $434.6567655$ & $386.3615693$ & $338.0663732$ & $289.771177$ & $241.4759808$\\ $5$ & $478.7815616$ & $430.9034055$ & $383.0252493$ & $335.1470931$ & $287.268937$ & $239.3907808$\\ $6$ & $474.647174$ & $427.1824566$ & $379.7177392$ & $332.2530218$ & $284.7883044$ & $237.323587$\\ $7$ & $470.5484877$ & $423.4936389$ & $376.4387902$ & $329.3839414$ & $282.3290926$ & $235.2742439$\\ $8$ & $466.4851945$ & $419.8366751$ & $373.1881556$ & $326.5396362$ & $279.8911167$ & $233.2425973$\\ $9$ & $462.4569888$ & $416.2112899$ & $369.965591$ & $323.7198922$ & $277.4741933$ & $231.2284944$\\ $10$ & $458.4635675$ & $412.6172108$ & $366.770854$ & $320.9244973$ & $275.0781405$ & $229.2317838$\\ $11$ & $454.5046304$ & $409.0541673$ & $363.6037043$ & $318.1532413$ & $272.7027782$ & $227.2523152$\\ $12$ & $450.5798795$ & $405.5218916$ & $360.4639036$ & $315.4059157$ & $270.3479277$ & $225.2899398$\\ $13$ & $446.6890198$ & $402.0201178$ & $357.3512159$ & $312.6823139$ & $268.0134119$ & $223.3445099$\\ $14$ & $442.8317586$ & $398.5485827$ & $354.2654069$ & $309.982231$ & $265.6990551$ & $221.4158793$\\ $15$ & $439.0078056$ & $395.1070251$ & $351.2062445$ & $307.3054639$ & $263.4046834$ & $219.5039028$\\ $16$ & $435.2168734$ & $391.6951861$ & $348.1734987$ & $304.6518114$ & $261.130124$ & $217.6084367$\\ $17$ & $431.4586767$ & $388.3128091$ & $345.1669414$ & $302.0210737$ & $258.875206$ & $215.7293384$\\ $18$ & $427.732933$ & $384.9596397$ & $342.1863464$ & $299.4130531$ & $256.6397598$ & $213.8664665$\\ $19$ & $424.0393618$ & $381.6354256$ & $339.2314895$ & $296.8275533$ & $254.4236171$ & $212.0196809$\\ $20$ & $420.3776855$ & $378.339917$ & $336.3021484$ & $294.2643799$ & $252.2266113$ & $210.1888428$\\ ........ & ........ & ........ & ......... & ......... & ........ & ......... \\ ........ & ........ & ........ & ......... & ......... & ........ & ......... \\ $3600$ & $1.37882E-11$ & $1.24094E-11$ & $1.10305E-11$ & $9.65173E-12$ & $8.27291E-12$ & $6.89409E-12$\\ \hline \end{tabular} \end{center} \end{table} \newpage \begin{thebibliography}{99} \bibitem{altrock2015mathematics} Altrock, P.M., Liu, L.L. and Michor, F., 2015. The mathematics of cancer: integrating quantitative models. Nature Reviews Cancer, 15(12), pp.730-745. \bibitem{alvarado2016metabolic} Alvarado, A. and Arce, I., 2016. Metabolic functions of the lung, disorders and associated pathologies. Journal of clinical medicine research, 8(10), p.689. \bibitem{Abazari} Abazari MA, Soltani M, Kashkooli FM. Targeted nano-sized drug delivery to heterogeneous solid tumor microvasculatures: Implications for immunoliposomes exhibiting bystander killing effect. Physics of Fluids. 2023;35(1):011905 \bibitem{azizi2024mathematical} Azizi, T., 2024. Mathematical modeling of cancer progression. AppliedMath, 4(3), pp.1065-1079. \bibitem{Idrees2023} Idrees, M., Alnahdi, A.S. and Jeelani, M.B., 2023. Mathematical Modeling of Breast Cancer Based on the Caputo–Fabrizio Fractal-Fractional Derivative. Fractal and Fractional, 7(11), p.805. \bibitem{Adam} Adam, J.A.; Bellomo, N. A Survey of Models for Tumor-Immune System Dynamics; Springer Science \& Business Media: Berlin/Heidelberg, Germany, 1997 \bibitem{Akman} Akman, T., Arendt, L.M., Geisler, J., Kristensen, V.N., Frigessi, A. and Köhn-Luque, A., 2024. Modeling of mouse experiments suggests that optimal anti-hormonal treatment for breast cancer is diet-dependent. Bulletin of Mathematical Biology, 86(4), p.42. \bibitem{bates2004tobacco} Bates, C. and Rowell, A., 2004. Tobacco Explained... The truth about the tobacco industry... in its own words. \bibitem{berzins1998role} Berzins, S.P., Boyd, R.L. and Miller, J.F., 1998. The role of the thymus and recent thymic migrants in the maintenance of the adult peripheral lymphocyte pool. The Journal of experimental medicine, 187(11), pp.1839-1848. \bibitem{boffetta2003contribution} Boffetta, P. and Nyberg, F., 2003. Contribution of environmental factors to cancer risk. British medical bulletin, 68(1), pp.71-94. \bibitem{burdett1996adjuvant} Burdett, S., Pignon, J.P., Tierney, J., Tribodet, H., Stewart, L., Le Pechoux, C., Aupérin, A., Le Chevalier, T., Stephens, R.J., Arriagada, R. and Higgins, J.P., 1996. Adjuvant chemotherapy for resected early‐stage non‐small cell lung cancer. Cochrane Database of Systematic Reviews, 2015(3). \bibitem{cooper2006small} Cooper, S. and Spiro, S.G., 2006. Small cell lung cancer: treatment review. Respirology, 11(3), pp.241-248. \bibitem{curry1980fiber} Curry, F.E. and Michel, C.C., 1980. A fiber matrix model of capillary permeability. Microvascular research, 20(1), pp.96-99. \bibitem{cummings2008administration} Cummings, J., 2008. The administration of medicines. Clinical Skills for Student Nurses: Theory, Practice and Reflection, p.284. \bibitem{Chen2013} Chen C, Baumann WT, Clarke R, Tyson JJ (2013) Modeling the estrogen receptor to growth factor receptor signaling switch in human breast cancer cells. FEBS Lett 587(20):3327–3334 \bibitem{Chen2014} Chen C, Baumann WT, Xing J, Xu L, Clarke R, Tyson JJ (2014) Mathematical models of the transitions between endocrine therapy responsive and resistant states in breast cancer. J R Soc Interface 11(96):20140206 \bibitem{van2016mitochondrial} Van den Bossche, J., Baardman, J., Otto, N.A., van der Velden, S., Neele, A.E., van den Berg, S.M., Luque-Martin, R., Chen, H.J., Boshuizen, M.C., Ahmed, M. and Hoeksema, M.A., 2016. Mitochondrial dysfunction prevents repolarization of inflammatory macrophages. Cell reports, 17(3), pp.684-696. \bibitem{devita2012two} DeVita Jr, V.T. and Rosenberg, S.A., 2012. Two hundred years of cancer research. New England Journal of Medicine, 366(23), pp.2207-2214. \bibitem{durovski2023insights} Durovski, D., Jankovic, M. and Prekovic, S., 2023. Insights into androgen receptor action in lung cancer. Endocrines, 4(2), pp.269-280. \bibitem{debela2021new} Debela, D.T., Muzazu, S.G., Heraro, K.D., Ndalama, M.T., Mesele, B.W., Haile, D.C., Kitui, S.K. and Manyazewal, T., 2021. New approaches and procedures for cancer treatment: Current perspectives. SAGE open medicine, 9, p.20503121211034366. \bibitem{dempsey2015distinct} Dempsey, L.A., 2015. Distinct tumor APCs. Nature Immunology, 16(1), pp.12-12. \bibitem{ding2025characterising} Ding, B., 2025. Characterising Tumour Volume and Growth Dynamics in Non-Small Cell Lung Cancer (Doctoral dissertation, UCL (University College London)). \bibitem{Doisneau2023} Doisneau-Sixou S, Sergio C, Carroll J, Hui R, Musgrove E, Sutherland R (2003) Estrogen and antiestrogen regulation of cell cycle progression in breast cancer cells. Endocr Relat Cancer 10(2):179–186 \bibitem{dwivedi2020perspective} Dwivedi, M., 2020. A perspective review on cancer–The deadliest disease. International Journal of Cancer, 1(1). \bibitem{eftimie2021mathematical} Eftimie, R. and Barelle, C., 2021. Mathematical investigation of innate immune responses to lung cancer: The role of macrophages with mixed phenotypes. Journal of Theoretical Biology, 524, p.110739. \bibitem{fialkow1976clonal} Fialkow, P.J., 1976. Clonal origin of human tumors. Biochimica et Biophysica Acta (BBA)-Reviews on Cancer, 458(3), pp.283-321. \bibitem{greim2018introduction} Greim, H. and Snyder, R., 2018. Introduction to the Discipline of Toxicology. Toxicology and Risk Assessment: A Comprehensive Introduction, pp.1-19. \bibitem{gomez2020heterogeneity} Gomez, H., 2020. How heterogeneity drives tumour growth: a computational study. Philosophical Transactions of the Royal Society A, 378(2171), p.20190244. \bibitem{harshe2023predicting} Harshe, I., Enderling, H. and Brady-Nicholls, R., 2023. Predicting patient-specific tumor dynamics: how many measurements are necessary?. Cancers, 15(5), p.1368. \bibitem{hethcote2000mathematics} Hethcote, H.W., 2000. The mathematics of infectious diseases. SIAM review, 42(4), pp.599-653. \bibitem{hsia2016lung} Hsia, C.C., Hyde, D.M. and Weibel, E.R., 2016. Lung structure and the intrinsic challenges of gas exchange. Comprehensive physiology, 6(2), pp.827-895. \bibitem{hsu2017estrogen} Hsu, L.H., Chu, N.M. and Kao, S.H., 2017. Estrogen, estrogen receptor and lung cancer. International journal of molecular sciences, 18(8), p.1713. \bibitem{holliday2008administration} Holliday, L. and Kierulff, C., 2008. Administration of medicines. Clinical Skills in Child Health Practice, p.161. \bibitem{Idrees2021} Idrees, M.; Sohail, A. Bio-algorithms for the modeling and simulation of cancer cells and the immune response. Bio-Algorithms Med-Syst. 2021, \bibitem{italiani2014monocytes} Italiani, P. and Boraschi, D., 2014. From monocytes to M1/M2 macrophages: phenotypical vs. functional differentiation. Frontiers in immunology, 5, p.514. %%\bibitem{jasinski2010influence} Jasiński, M., 2010. Influence of emissivity changes on the blood flow rate determined on the basis of heat balance equation. Scientific Research of the Institute of Mathematics and Computer Science, 9(1), pp.37-44. \bibitem{jia2020study} Jia, B., Zhang, X., Mo, Y., Chen, B., Long, H., Rong, T. and Su, X., 2020. The study of tumor volume as a prognostic factor in T staging system for non-small cell lung cancer: an exploratory study. Technology in Cancer Research \& Treatment, 19, p.1533033820980106. \bibitem{juan2021chemistry} Juan, C.A., Pérez de la Lastra, J.M., Plou, F.J. and Pérez-Lebeña, E., 2021. The chemistry of reactive oxygen species (ROS) revisited: outlining their role in biological macromolecules (DNA, lipids and proteins) and induced pathologies. International journal of molecular sciences, 22(9), p.4642. \bibitem{Tufail} Tufail, M., Cui, J. and Wu, C., 2022. Breast cancer: molecular mechanisms of underlying resistance and therapeutic approaches. American journal of cancer research, 12(7), p.2920. \bibitem{Carrillo} Ku-Carrillo RA, Delgadillo SE, Chen-Charpentier B (2016) A mathematical model for the effect of obesity on cancer growth and on the immune system response. Appl Math Model 40(7–8):4908–4920 %\bibitem{sengupta2021principles} SenGupta, S., Parent, C.A. and Bear, J.E., 2021. The principles of directed cell migration. Nature Reviews Molecular Cell Biology, 22(8), pp.529-547. \bibitem{Simpson} Simpson, E.R. and Dowsett, M., 2002. Aromatase and its inhibitors: significance for breast cancer therapy. Recent progress in hormone research, 57, pp.317-338. \bibitem{siegfried2009estrogen} Siegfried, J.M., Hershberger, P.A. and Stabile, L.P., 2009, December. Estrogen receptor signaling in lung cancer. In Seminars in oncology (Vol. 36, No. 6, pp. 524-531). WB Saunders. \bibitem{spiro2010lung} Spiro, S.G., Tanner, N.T., Silvestri, G.A., Janes, S.M., Lim, E., Vansteenkiste, J.F. and Pirker, R., 2010. Lung cancer: progress in diagnosis, staging and therapy. Respirology, 15(1), pp.44-50. \bibitem{mcadam2016influence} McAdam, K., Eldridge, A., Fearon, I.M., Liu, C., Manson, A., Murphy, J. and Porter, A., 2016. Influence of cigarette circumference on smoke chemistry, biological activity, and smoking behaviour. Regulatory Toxicology and Pharmacology, 82, pp.111-126. \bibitem{martin2022role} Martin-Perez, M., Urdiroz-Urricelqui, U., Bigas, C. and Benitah, S.A., 2022. The role of lipids in cancer progression and metastasis. Cell metabolism, 34(11), pp.1675-1699. %\bibitem{Martin} Martin, H.L., Smith, L. and Tomlinson, D.C., 2014. Multidrug-resistant breast cancer: current perspectives. Breast Cancer: targets and therapy, pp.1-13. \bibitem{manukonda2025dose} Manukonda, M.R.M. and Scholar, P.G., 2025. DOSE-RESPONSE RELATIONSHIPS. Principles of Medical Toxicology, p.36. \bibitem{maitra2021targeting} Maitra, R., Malik, P. and Mukherjee, T.K., 2021. Targeting estrogens and various estrogen-related receptors against non-small cell lung cancers: a perspective. Cancers, 14(1), p.80. \bibitem{mann2009principles} Mann, U., 2009. Principles of chemical reactor analysis and design: new tools for industrial chemical reactor operations. John Wiley \& Sons. \bibitem{mc2012historical} Mc Laughlin, J., 2012. An historical overview of radon and its progeny: applications and health effects. Radiation protection dosimetry, 152(1-3), pp.2-8. \bibitem{mital2007thermal} Mital, M. and Scott, E.P., 2007. Thermal detection of embedded tumors using infrared imaging. Journal of biomechanical engineering, 129(1), pp.33-39. \bibitem{nair2012isolation} Nair, S., Archer, G.E. and Tedder, T.F., 2012. Isolation and generation of human dendritic cells. Current protocols in immunology, 99(1), pp.7-32. \bibitem{Normanno} Normanno N, Di Maio M, De Maio E, De Luca A, De Matteis A, Giordano A, Perrone F (2005) Mechanisms of endocrine resistance and novel therapeutic strategies in breast cancer. Endocr Relat Cancer 12(4):721–747 \bibitem{ni1997role} Ni, K. and O'neill, H.C., 1997. The role of dendritic cells in T cell activation. Immunology and cell biology, 75(3), pp.223-230. %\bibitem{ndreko2025modeling} Ndreko, E., Victor, S., Duhé, J.F. and Melchior, P., 2025. Modeling of bio-heat transfers in lungs with fractional models. Annual Reviews in Control, 60, p.101010. \bibitem{ozkose2021fractional} Özköse, F., Yılmaz, S., Yavuz, M., Öztürk, İ., Şenel, M.T., Bağcı, B.Ş., Doğan, M. and Önal, Ö., 2021. A fractional modeling of tumor–immune system interaction related to lung cancer with real data. The European Physical Journal Plus, 137(1), p.40. \bibitem{pandi2016brief} Pandi, A., Mamo, G., Getachew, D., Lemessa, F., Kalappan, V.M. and Dhiravidamani, S., 2016. A brief review on lung cancer. Int. J. Pharma Res. Health Sci, 4, pp.907-914. \bibitem{patel2017fate} Patel, A.A., Zhang, Y., Fullerton, J.N., Boelen, L., Rongvaux, A., Maini, A.A., Bigley, V., Flavell, R.A., Gilroy, D.W., Asquith, B. and Macallan, D., 2017. The fate and lifespan of human monocyte subsets in steady state and systemic inflammation. Journal of Experimental Medicine, 214(7), pp.1913-1923. %\bibitem{palucka2012cancer} Palucka, K. and Banchereau, J., 2012. Cancer immunotherapy via dendritic cells. Nature Reviews Cancer, 12(4), pp.265-277. \bibitem{peng2016computational} Peng, H., Tan, H., Zhao, W., Jin, G., Sharma, S., Xing, F., Watabe, K. and Zhou, X., 2016. Computational systems biology in cancer brain metastasis. Frontiers in bioscience (Scholar edition), 8(1), p.169. \bibitem{Pontryagin} Pontryagin LS, Boltyanskii VT, Gamkrelidze RV, Mishcheuko EF. The mathematical theory of optimal processes, Wiley, New Jersey, 1962. \bibitem{proctor2012history} Proctor, R.N., 2012. The history of the discovery of the cigarette–lung cancer link: evidentiary traditions, corporate denial, global toll. Tobacco control, 21(2), pp.87-91. \bibitem{lai2018modeling}Lai, X., Stiff, A., Duggan, M., Wesolowski, R., Carson III, W.E. and Friedman, A., 2018. Modeling combination therapy for breast cancer with BET and immune checkpoint inhibitors. Proceedings of the National Academy of Sciences, 115(21), pp.5534-5539. \bibitem{lu2018haart}Lu, D.Y., Wu, H.Y., Yarla, N.S., Xu, B., Ding, J. and Lu, T.R., 2018. HAART in HIV/AIDS treatments: future trends. Infectious Disorders-Drug TargetsDisorders), 18(1), pp.15-22. \bibitem{Luque} Luque-Bolivar, A., Pérez-Mora, E., Villegas, V.E. and Rondón-Lagos, M., 2020. Resistance and overcoming resistance in breast cancer. Breast Cancer: Targets and Therapy, pp.211-229. \bibitem{rosenberg2007erwin} Rosenberg, C.E., 2007. Erwin H. Ackerknecht, social medicine, and the history of medicine. Bulletin of the History of Medicine, 81(3), pp.511-532. \bibitem{roszkowska2024multilevel} Roszkowska, M., 2024. Multilevel mechanisms of cancer drug resistance. International journal of molecular sciences, 25(22), p.12402. \bibitem{saini2020cancer} Saini, A., Kumar, M., Bhatt, S., Saini, V. and Malik, A., 2020. Cancer causes and treatments. Int J Pharm Sci Res, 11(7), pp.3121-3134. \bibitem{singh2022mathematical} Singh, K., 2022. Mathematical Modelling of Metastatic Processes on a Network (Master's thesis, University of the Witwatersrand, Johannesburg (South Africa)). \bibitem{simoes2015metabolic} Simões, R.V., Serganova, I.S., Kruchevsky, N., Leftin, A., Shestov, A.A., Thaler, H.T., Sukenick, G., Locasale, J.W., Blasberg, R.G., Koutcher, J.A. and Ackerstaff, E., 2015. Metabolic plasticity of metastatic breast cancer cells: adaptation to changes in the microenvironment. Neoplasia, 17(8), pp.671-684. \bibitem{sabir2025mathematical} Sabir, S., León-Triana, O., Serrano, S., Barrio, R. and Pérez-García, V.M., 2025. Mathematical model of CAR T-cell therapy for a B-cell lymphoma lymph node. Bulletin of Mathematical Biology, 87(3), pp.1-33. \bibitem{spiro2005one} Spiro, S.G. and Silvestri, G.A., 2005. One hundred years of lung cancer. American journal of respiratory and critical care medicine, 172(5), pp.523-529. \bibitem{Sun} Sun, X., Bao, J. and Shao, Y., 2016. Mathematical modeling of therapy-induced cancer drug resistance: connecting cancer mechanisms to population survival rates. Scientific reports, 6(1), p.22498. %\bibitem{shadden2015lagrangian} Shadden, S.C. and Arzani, A., 2015. Lagrangian postprocessing of computational hemodynamics. Annals of biomedical engineering, 43(1), pp.41-58. \bibitem{tallarida2012dose} Tallarida, R.J. and Jacob, L.S., 2012. The dose—Response relation in pharmacology. Springer Science \& Business Media. \bibitem{tao2019epidemiology} Tao, M.H., 2019. Epidemiology of lung cancer. Lung Cancer and Imaging, pp.4-1. \bibitem{tarin2011cell} Tarin, D., 2011, April. Cell and tissue interactions in carcinogenesis and metastasis and their clinical significance. In Seminars in cancer biology (Vol. 21, No. 2, pp. 72-82). Academic Press. \bibitem{talkington2015estimating} Talkington, A. and Durrett, R., 2015. Estimating tumor growth rates in vivo. Bulletin of mathematical biology, 77(10), pp.1934-1954. \bibitem{tripathy2020cancer} Tripathy, A. and Pradhan, R.K., 2020. The Cancer Principle-III: Evolution by Cancer. Indian Journal of Public Health Research \& Development, 11(6). \bibitem{timmermann2014lung} Timmermann, C., 2014. Lung cancer and consumption in the nineteenth century: bodies, tissues, cells and the making of a rare disease. In A History of Lung Cancer: The Recalcitrant Disease (pp. 11-33). London: Palgrave Macmillan UK. \bibitem{uthamacumaran2022review} Uthamacumaran, A. and Zenil, H., 2022. A review of mathematical and computational methods in cancer dynamics. Frontiers in oncology, 12, p.850731. \bibitem{van2016mitochondrial} Van den Bossche, J., Baardman, J., Otto, N.A., van der Velden, S., Neele, A.E., van den Berg, S.M., Luque-Martin, R., Chen, H.J., Boshuizen, M.C., Ahmed, M. and Hoeksema, M.A., 2016. Mitochondrial dysfunction prevents repolarization of inflammatory macrophages. Cell reports, 17(3), pp.684-696. \bibitem{villa2025reducing} Villa, C., Maini, P.K., Browning, A.P., Jenner, A.L., Hamis, S. and Cassidy, T., 2025. Reducing phenotype-structured partial differential equations models of cancer evolution to systems of ordinary differential equations: a generalised moment dynamics approach. Journal of Mathematical Biology, 91(2), p.22. \bibitem{venter2004century} Venter, C. and Cohen, D., 2004. The century of biology. New Persp. Q., 21, p.73. \bibitem{vineis2014global} Vineis, P. and Wild, C.P., 2014. Global cancer patterns: causes and prevention. The Lancet, 383(9916), pp.549-557. \bibitem{wang2025optimal} Wang, D. and Lei, J., 2025. Optimal adaptive therapeutic schedules for metastatic castrate-resistant prostate cancer based on bilevel optimization problem. Journal of Mathematical Biology, 90(6), p.60. %\bibitem{welf2011signaling} Welf, E.S. and Haugh, J.M., 2011. Signaling pathways that control cell migration: models and analysis. Wiley Interdisciplinary Reviews: Systems Biology and Medicine, 3(2), pp.231-240. \bibitem{xiong2025mathematical} Xiong, Z., Xia, Y., Xue, L. and Lei, J., 2025. Mathematical Modelling and Optimization of Medication Regimens for Combination Immunotherapy of Breast Cancer: Z. Xiong et al. Bulletin of Mathematical Biology, 87(7), p.88. \bibitem{yu2016occurrence} Yu, Y., Cui, Y., Niedernhofer, L.J. and Wang, Y., 2016. Occurrence, biological consequences, and human health relevance of oxidative stress-induced DNA damage. Chemical research in toxicology, 29(12), pp.2008-2039. \bibitem{zhang2022evolution} Zhang, J., Cunningham, J., Brown, J. and Gatenby, R., 2022. Evolution-based mathematical models significantly prolong response to abiraterone in metastatic castrate-resistant prostate cancer and identify strategies to further improve outcomes. Elife, 11, p.e76284. %\bibitem{Lorains} % %\bibitem{Bu} \end{thebibliography} \end{document}