Summary five major difficulties of modern semiconductor manufacturing.
\r\n\tThis book will provide an overview of current knowledge regarding the role of AMPK and mechanisms through which it contributes to energy homeostasis in mammalian physiology and pathology. Topics will broadly include, but not be limited to:
\r\n\tMechanisms of AMPK regulation,
\r\n\tLipid and carbohydrate metabolism,
\r\n\tCellular energy sensing and homeostasis,
\r\n\tAdaptation to exercise and fasting,
\r\n\tMitochondrial physiology, pathology, and dynamics,
\r\n\tLifespan and ageing,
\r\n\tTissue injury and organ failure,
The semiconductor plant investment exceeds 3 billion US dollars, and the basic operating cost per day exceeds 6 million US dollars, although this process is very important, especially in commercial key sizes below 12 nm (Figure 1) [1, 2]. However, modern semiconductor manufacturing has five major difficulties, which are summarized in Table 1. Semiconductor refers to the material whose conductivity is between conductor and insulator at normal temperature. Semiconductors are used in integrated circuits, consumer electronics, communication systems, photovoltaic power generation, lighting applications, high-power power conversion and other fields. Such as diode is a device made of semiconductor. Whether from the perspective of technology or economic development, the importance of semiconductors is very huge [3, 4]. 5G and artificial intelligence technologies are the main application areas of advanced technology. For example, in 2019, many mobile phone manufacturers launched 5G models, and most of these models used 12.01 ∼ 7.01 advanced technology baseband chips, such as Qualcomm Snapdragon X50, MediaTek Helio M70, Intel XMM8000 series, Samsung Exynos 5000 series, Hisilicon Balong 5000 series, etc. This creates a potential demand for semiconductor equipment. 5G and artificial intelligence will not only bring the semiconductor equipment market back to a short term in 2020, but also support the development of the semiconductor equipment industry in the long term. The research and development organization predicts the technology trend in 2020, and the research and development of Toyo, a subsidiary of Jibang, predicts that the semiconductor industry will gradually come out of the bottom in 2020 with the continuous increase in demand for 5G, AI, automotive and other emerging applications. Among them, 5G transformative technologies are the most critical. The Industrial Intelligence Institute (MIC) of the China Resources Planning Association pointed out that this year, about 56 telecom operators in 32 countries have announced the deployment of 5G networks, of which 39 telecom operators have officially opened 5G services. It is estimated that by 2020, there will be 170 telecom providers worldwide providing 5G commercial services [5, 6].
|Integration is getting higher and higher||The number of transistors integrated on a chip is increasing, from the 1960s to the present, from one transistor to more than 10 billion.|
|The higher the accuracy requirements, the more difficult the process.||The key size was reduced from 1um in 1988 to 5 nm in 2020, a reduction of 99.5%. From this perspective, the difficulty of integrated circuit manufacturing is gradually increasing, and the acceleration of the difficulty is also increasing.|
|Difficult to break through some specific single-point technology||The process technology that constitutes the smallest unit of the semiconductor manufacturing process is the single-point technology, especially the photolithography process. The manufacturing process of complex circuits exceeds 500 processes. These processes are carried out under precision instruments, which are not clear to human eyes.|
|Need to integrate multiple technologies||The difficulty of integrating technology lies in how to complete a technological process with low cost, satisfactory specifications and complete operation from an unlimited combination of component technologies in a short time. This includes operators, process equipment, raw materials, process parameters and methods, FAB environmental conditions and so on.|
|Mass production technology||The process flow constructed by the R & D center through integrated technology is transferred to the mass production plant. Strict exact copying is basically impossible. Even if the equipment in the development center and the mass production plant are the same, the same result may not be obtained under the same process conditions. This is because even with the same equipment, there will be a slight performance difference between the two machines. Therefore, with the continuous improvement of the degree of semiconductor precision, the problem of machine difference is becoming more and more obvious.|
Combining the foregoing, in order to achieve high accuracy and high throughput, modern FAB uses a large number of high-energy processes, such as plasma, CVD, and ion implantation. This furnace is one of the important tools for semiconductor manufacturing (Figure 2). Due to the high energy, the physicochemical changes in each reactor are very complex, and it is often impossible to determine the type and concentration of the by-products produced, and the type and concentration of these by-products often change randomly. These by-products usually cause the following effects on FAB, including: (1) Incompatibility between by-products may increase the toxicity or explosiveness of the gas in the pipeline; (2) By-products may corrode the exhaust pipe or make it brittle (3) If the type and concentration of by-products cannot be determined, it is impossible to select the appropriate exhaust gas treatment equipment; (4) The currently used processing equipment may be damaged, which may affect the processing efficiency [7, 8, 9]. As of 2020, there are already 37 semiconductor wafers FABs in Taiwan. Taking the 12-inch FAB as an example, there are about 420 various types of main process machines in a manufacturing plant with a monthly capacity of 50,000 wafers, of which about 63 are thin film process machines, photolithography process machine also has about 55, summary of Taiwan IC manufacturing FABs show in Table 2 .
|Trade names||5″ FAB||6″ FAB||8″ FAB||12″ FAB|
|Micron||—||—||—||3 (Taichung * 1 + Hua Ya * 2)|
|MXIC 2||—||1||1||1 (Hsinchu Science Park)|
|NanYa||—||—||1||1 (Taishan Nanlin Park)|
|PSC||—||—||—||3 (Hsinchu Science Park)|
|TSMC||—||1||7||5 (Hsinchu Science Park *3 + Tainan Science Park*2)|
|UMC||—||1||6||3 (Hsinchu Science Park)|
|Winbond||—||—||—||2 (Taichung Science Park +Tainan Science Park)|
Treat 100,000 pieces of factory space as a large factory, that is, all machines with similar functions are only divided into a group. The layout of the machine group is planned using a typical “non-shaped” pattern (please refer Figure 2(b)). A FAB is divided into several rectangular blocks, and each block is called a processing area (bay). Machines in the same machine group or machines with similar functions are placed on the same bay as the principle. The bay is located on the two sides of the FAB. In the center is the among-bay goods material handling system, which is responsible for the transportation between bay and bay. Within each bay, there is also a within-bay goods material handling system, which is responsible for the machine and transport with the machine .
However, due to the FAB12″ 7 nm stove based on the above production management requirements, the FTIR system was installed in this study, which includes (1) confirming the characteristics of harmful process exhaust gas; (2) evaluating the processing efficiency of various process equipment process exhaust; (3) Conduct a hazard exposure assessment survey during machine maintenance and repair; (4) Confirm the concentration and source of harmful gases and particles in the clean room operating environment; (5) Identify harmful substances in the pipeline. It is intelligent with the hardware system, and the IoT module is added to the original module. In this study, various process parameters and information required by FAB are continuously obtained in the 12″ stove [12, 13].
Under this premise, in order to make the above software and hardware system intelligent, add the IoT module to the original module, so that we can continuously obtain various process parameters and information required by FAB in the 12″ furnace process tool. Through 24 hr of continuous Processing and thousands of processing machines, we will obtain a large amount of data to confirm the above production requirements, so that we can effectively master the FAB characteristics, improve production efficiency, improve product yield and establish a safe and healthy product line and employee working environment It is an important contribution of this research [14, 15, 16].
The core of the FTIR is the Michealson interferometer. Its principle is that the two infrared beams after the infrared light source is split by the beam splitter are respectively directed to the fixed mirror and the moving mirror, and then combined into A single infrared ray, due to the difference in optical path formed by the moving mirror, makes the final combined infrared ray form an infrared beam of different energy due to destructive and constructive interference. It has a fast analysis speed, is not destructive to the sample, and can analyze solid liquid and gas samples, making it gradually become an indispensable qualitative tool for material analysis. In certain circumstances, it can even achieve the ability to quickly screen quantitative .
The instrument used in this study is aspirated FTIR. The pumped FTIR uses a pump to introduce the gas to be tested into the FTIR detection chamber for immediate analysis. The measurement method is shown in Figure 3. The main components of the pumped FTIR include infrared sources, interferometers, beam splitters, fixed mirrors, moving mirrors and gas chambers, detectors and electronic modules, etc. In addition, there must be a sampling tube and pump and other gas samples into the closed cavity In addition to the computer used for data collection and data analysis and appropriate software, this study also added an IoT module to the existing FTIR, allowing the FTIR to transmit and calculate. With the cloud, pumped FTIR’s the instrument configuration is shown in Figure 4 [18, 19]. For IoT part, this study starts with the sensor, imports the sensing signal into the electronic module, and exports the signal to the cloud system through the WiFi module, in this way, this study obtain FTIR sensing data 24 hr, and this system is set to obtain data once per second.
The basic design of the infrared spectrometer is to emit a beam of light to the measurement area and measure the amount of intensity change after the beam passes through the gas to be tested. Since each gas molecule has its specific infrared light absorption coefficient, when a light beam passes through the measurement region, a specific gas molecule absorbs light of a specific wavelength, so that the intensity of the light beam in this wavelength band is weakened, and the ratio of light intensity before and after absorption is The concentration of the gas is directly related. The absorption band and intensity of the gas sample can be measured to know the composition and concentration contained in the gas. For a maximum path difference d adjacent wavelengths λ1 and λ2 will have n and (n + 1) cycles respectively in the interferogram. The corresponding frequencies are ν1 and ν2, and the membership function in the following Eqs. (1)∼(5) [20, 21]:
FTIR mainly emits a beam of light to the measurement area and measures the intensity change of the beam after passing the gas to be measured. Since each gas molecule has its specific infrared light absorption coefficient, when the light beam passes through the measurement area, the specific gas molecule will absorb light of a specific wavelength, so that the intensity of the light beam in this band is reduced, and the ratio of the light intensity before and after absorption The concentration of the gas is directly related, and the absorption band and intensity of the gas sample can be calculated to know the composition and concentration of the gas. Figure 5 shows the absorption spectrum of several gas molecules in the infrared range.
This research is based on the fact that FTIR operates continuously for 24 hr and captures signals every second to obtain data throughout the year as a basis for big data. The main research significance is not to grasp the huge data information, but to professionally process these meaningful data, so that FAB factory managers can know whether the exhaust emissions meet the alert or warning potential under the sensing trend. Most of the research uses the big data platform for deployment, debugging and maintenance. This study is connected to the Chief Cloud eXchange (CCX) cloud database platform and the Spark big data platform. At the same time, these two platforms are also more suitable for primary systems to implement systems. Figure 6 is the research cloud computing system and database architecture.
For the big data of the chip process exhaust obtained by RFID, full consideration should be given to (1) big data life cycle, (2) big data technology ecology, (3) big data acquisition and preprocessing, (4) big data storage and management, (5) Big data computing model and system.
The big data analysis method used in this study is described as follows. U is defined as the non-empty initial universe of the object. Then define E as a set of parameters related to the object in U. Let P (U) be the power set of U, and A ⊂ E. A pair (F, A) is called a soft set on U, where F is the mapping given by a F: A → (U). In other words, the soft set on U is a parameterized family of Universe U subsets.
In addition, if the universe set U is a non-empty finite set, and σ is the equivalent relationship on U. Then (U, σ) is called approximate space. If X is a subset of U, then X can be written as a union of equivalent classes of U or not. If X can be written as a union of equivalent classes of U, then X is definable, otherwise it is undefinable. If X is undefinable, it can be approximated as two definable subsets, called the upper and lower approximation of X, as shown below Eq. (6) [22, 23].
If X contains remainder, then
Where the number of groups will be added to 1.
Algorithm (1): The most optimized attribute set searching algorithm.
Input: Optimized reduct sets, R1 until Rn Output: The most optimal reduct set.
if Reduct set R has more than one value then.
Select the highest number of attribute values, HR if HR does not have the same number with attribute value AND HR has more than one value then.
Select the first reduction set, FR of attribute values.
Proceed to the next process
Proceed to the next process
Algorithm (2): Soft set parameter reduction algorithm.
In tabular representation, let (F, P) represent the soft set. If Q is the reduction of P, the soft set reduction set is defined as (F, Q) of the soft set (F, P) where P ⊂ E.
Input: A soft set (F, E), set P.
Output: Optimal decision.
Input the set P of choice parameters.
Find all reducts of (F, P).
Select one reduct set (F, Q) of (F, P).
Find weighted table of soft set (F, Q) according to the decided weights.
Find k, for which ck = max ci.
hk is the optimal choice of value for the selected object. If k has more than one value, any one of the benefits could be chosen.
ci is the choice of value of an object hi where and hi,j is the entries in the table of the reduct soft set.
Algorithm (3): Rough set parameter reduction algorithm.
Input: An information system S = (U, A, V, f).
U is a finite nonempty set object.
A is a finite nonempty set of attributes.
V is a nonempty set of values.
f is an information function that maps an.
object in U to exactly one value in V.
Output: Simplified reduct sets.
Input the information Table S.
Discretization of data.
Forming up the n × n discernibility matrix. The elements of S table is defined as d(x, y) = a ∈ A | f (x, a) ̸= f (y, a), d(x, y) is an attribute.
set distinguishing x and y. For each attribute a ∈ A, if d(x, y) = a1, a2, …, ak ̸= ∅.
Formulate the Boolean function a1 ∨ a2… ∨ ak or discernibility function which represented by Σ d(x, y) as indicated: F (A) = Π (x,y)∈U × U) Σ d(x, y).
If d(x, y) = ∅, constant 1 will be assigned to the Boolean function.
Execute the attribute reduction process based on the simplified Boolean function.
New optimized reduct sets are generated.
In the reactor, chemical reaction is used to form the reactant (usually a gas) into a solid product, and a thin film is deposited on the surface of the wafer. This process is called CVD (Chemical Vapor Deposition). This process has (1) good step coverage, (2) energy with high aspect ratio gap filling, (3) good thickness uniformity, (4) high pure and dense, (5) when ratio can be controlled, (6) low stress for high film quality, (7) good electrical properties, (8) base plate and excellent film adhesion characteristics. Figure 7 shows the system architecture of CVD in FAB.
The precipitation of products during the CVD process can be divided into the following steps: (1) the source gas diffuses to the substrate surface, (2) the substrate adsorbs the source gas, (3) the substances adsorbed on the substrate react chemically on its surface, (4) The precipitated material diffuses on the surface of the substrate, (5) The reaction product is separated from the gas-phase reactant, (6) The precipitated non-volatile material is deposited on the substrate surface by diffusion and the like. The composition, structure and performance of the products obtained in this chemical reaction can be controlled by changing the parameters of the reaction. The reaction parameters mainly include the type of gas, the gas reaction concentration, the delivery method of the reactant, the gas flow rate, the total gas pressure and the Area pressure, heating method, substrate material, substrate surface state, substrate reaction temperature, temperature distribution and gradient, etc.
When entering the 7-nanometer process, the channel material of the semiconductor PN junction must also be changed. Since the electron mobility of silicon is 1500 cm2 / Vs, and germanium can reach 3900 cm2 / Vs, and the implementation voltage of silicon devices is 0.75 ∼ 0.8 V, while the germanium devices are only 0.5 V, so germanium was Considered to be the preferred material for MOSFET transistors, the first 7-nanometer wafer in IBM Lab used Ge-Si material. IMEC researched new germanium-doped materials and screened two channel materials that can be used for 7 nm: one is PFET composed of 80% germanium and the other is 25 ∼ 50% mixed germanium FET Or 0 ∼ 25% NFET mixed with germanium (Figure 8) [24, 25].
The silicon dioxide in the device, in this study, uses TEOS as the raw material, Tetraethyl Orthosilicate, the chemical formula is Si (OC2H5)4. High boiling point (about 169° C under normal pressure), store and use in liquid form. TEOS is liquid at room temperature and normal pressure. In order to increase the use of CVD process and the stability of the process, the TEOS container (about 40 ∼ 70°C) is heated during use to increase its saturated vapor pressure In the gaseous use of TEOS in the deposition reaction of CVD. The process parameters and characteristics of TEOS are summarized in Table 3.
|Deposition temperature (°C)||650 ∼ 750||300 ∼ 400||350 ∼ 450|
|Operating pressure (Torr)||1 ∼ 10||0.1 ∼ 5||500 ∼ 700|
|Refractive index||1.43 ∼ 1.46||1.47 ∼ 1.5||1.45|
|Dielectric constant||4.0||4.1 ∼ 4.9||4.4|
|BOE(100:1) Etching rate ( )||30||400||1200|
|Stress value (dyne/cm2)||1 ∼ 3 × 109||−(1 ∼ 5 × 109)||108 ∼ 3 × 109|
According to the “Semiconductor Manufacturing Air Pollution Control and Emission Standards” announced by the Environmental Protection Agency of the Taiwan Government, air pollutants produced in the process should be discharged after being purified by the appropriate system, where the efficiency of the system or the total emissions of the entire factory should be Meet the standards listed in the Table 4.
|Air pollutants||Equipment efficiency standard||Total control standards|
|Volatile Organic Compounds||>90%||<0.6 kg/hr (calculated based on methane)|
|Nitric acid, hydrochloric acid, phosphoric acid and hydrofluoric acid||>95%||<0.6 kg/hr|
|Sulfuric acid droplets||>95%||<0.1 kg/hr|
In order to achieve high precision and high output, modern high-energy processes such as plasma are important tools for semiconductor manufacturing. The physicochemical changes that occur in each reactor due to high energy are quite complicated, and the type and concentration of by-products cannot often be determined. These by-products usually have the following effects on the plant. (1) Incompatibility between by-products may increase the toxicity or explosiveness of the gas in the pipeline; (2) By-products may cause corrosion or embrittlement of exhaust pipe materials; (3) If the type and concentration of by-products cannot be determined, appropriate exhaust gas treatment equipment may be selected; (4) damage may be caused to the currently used treatment equipment, which may affect treatment efficiency. Based on the aforementioned production management requirements, the 12″ furnace FTIR system installed by FAB includes (1) confirmation of the characteristics of hazardous process exhaust gas; (2) evaluation of the processing efficiency of various process exhaust gas treatment equipment; (3) investigation of hazardous sources. Condition assessment during machine maintenance and repair; (4) Confirm the concentration and source of harmful gases and particulate matter in the clean room operating environment [26, 27, 28].
In this study, open-path FTIR was used to monitor the air quality of clean room to ensure the air quality of the working environment and the health of employees. The principle of measurement is the same as the principle of Extractive FTIR, but the closed cavity is changed to open type and integrated IoT mechanism to connect to the cloud, which is suitable for a variety of gaseous pollutants (including organic gases and inorganic gaseous pollutants) in the atmosphere are monitored (Figure 9). Figure 10 shows the situation where this study was set on the site to set the exhaust line of the furnace control process.
This study set up two measuring points in the 12″ factory of Hsinchu Science Park in Taiwan, as shown in Figure 11. Table 5 shows the process parameters of the on-site process tools during our experiment, and Table 6 shows the processing parameters of the on-site machine exhaust gas treatment equipment. Among them, Inlet flow rate (Qi) estimated from the TEOS injection = 89 LPM, Initial outlet flow rate (Qo) estimated from the TEOS injection = 292 LPM. Therefore, dilution ratio = Qo/Qi = 292/89 = 3.3.
|Process tools||BPSG of A point||BPSG of b point|
|TMB flow||30 sccm||30 sccm|
|TEOS flow||300 sccm||300 sccm|
|PH3 flow||0.77 slm||0.8 slm|
|Chamber pressure||1.1 torr||0.8 torr|
|concentration compound||Inlet||Outlet||Efficiency (%)|
|Max (ppm)||Average (ppm)||Max (ppm)||Average (ppm)|
From the measurement data of this study, it can be found (Figure 12) that the main reactant of the thin film process is TEOS for BPSG, so almost all reactions are carried out in the reaction chamber, or become a composite, which does not exist in the main exhaust gas pipeline, the beginning of the reaction concentration between 0.08 ∼ 0.1 ppm only. C2H4 is mainly used for cleaning the reaction chamber, concentration between 0.15 ∼ 0.25 ppm. Therefore, when the main process is carried out, the high concentration of the input will be cleaned, so the high concentration state can be seen in the main exhaust pipe. On the other hand, the CO is active because it is in the process of production, so it is difficult to find the concentration of the main exhaust pipe. Figures 13 and 14 are the results of the study after optimizing the concentration values measured at the day of point A.
Figure 15 shows secondary main exhaust pipe concentration trend. Due to the proximity of the process chamber, the concentrations are clearly detected, especially CO is more obvious, concentration between 1 ∼ 1.7 ppm, and TEOS is liquid because it is normal and concentration between 0.4 ∼ 0.6 ppm, so although the concentration near the reaction chamber is high, condensation occurs when entering the low temperature zone, so only this The section pipeline is measured, and the main exhaust pipe is not obvious. The C2H4 concentration is between 1.5 ∼ 2.0 ppm.
According to the OHSA regulations, this study is set in the cloud database for big data analysis and decision making, when the upper limit of TEOS, C2H4, CO are 0.6, 2.0, 1.7 ppm; the lower limit of TEOS, C2H4, CO is 0.4, 1.5, 1 ppm. The application architecture of this study can be extended to other semiconductor processes, so that IoT integration and big data operations can be performed for all processes, this is an important step in promoting FAB intelligent production and an important contribution of this study.
In order to achieve high precision and yield, modern FABs use a large number of high-energy processes such as plasma, CVD and ion implantation, the furnace is one of the important tools of semiconductor manufacturing. The FAB installed FTIR system due to the 12″ furnace tools based on the aforementioned production management requirements. The principle of measurement is the same as the principle of Extractive FTIR, but the closed cavity is changed to open type and integrated IoT mechanism to connect to the cloud, which is suitable for a variety of gaseous pollutants (including organic gases and inorganic gaseous pollutants) in the atmosphere are monitored in this study. This study set up two measuring points of furnace process tools in the 12″ factory of Hsinchu Science Park in Taiwan. This study obtained FTIR measurements, and according to the OHSA regulations, this study is set in the cloud database for big data analysis and decision making, when the upper limit of TEOS, C2H4, CO are 0.6, 2.0, 1.7 ppm; the lower limit of TEOS, C2H4, CO is 0.4 , 1.5, 1 ppm. The application architecture of this study can be extended to other semiconductor processes, so that IoT integration and big data operations can be performed for all processes, this is an important step in promoting FAB intelligent production and an important contribution of this study.
The normal operation of the pituitary-thyroid axis depends on the levels of thyroid stimulating hormone (TSH) and thyroid hormones, triiodothyronine (T3) and thyroxine (T4) . Serum TSH is produced and released by the pituitary gland in response to the low levels of free thyroid hormones in the serum. Circulating TSH in turn stimulates the thyroid to produce and secrete T3 and T4 into the serum. When thyroid hormones reach their highest levels, it inhibits the production of TSH, which describes the normal day to day operation of the pituitary-thyroid axis. The axis is commonly referred as an important negative feedback loop in the endocrine system. The normal function of this loop is essential for the body’s metabolic rate—it affects how quickly the cells in our body use the energy stored within itself [2, 3]. For the purpose of this modeling work, we use only serum free T4 for the levels of thyroid hormones as frequently asked to measure in the thyroid clinics for all patients with Hashimoto’s autoimmune thyroiditis (see \nFigure 1\n).\n\n
The normal operation of the axis can be tested and verified clinically via one measurement of TSH and free T4 from the blood serum. Laboratory tests are important in diagnosing conditions of the thyroid gland . The result of the blood test is determined approximately based on the normal reference range of TSH and free T4 is (0.4–4) mU/L and (7–18) pg/mL, respectively as recommended by the American Thyroid Association . As the physiology of the axis is governed by TSH and free T4, the clinical state of the axis has been defined based on these values . Suppose TSH and free T4 levels falls within the normal reference range; the state of the axis is said to be clinically normal. Suppose TSH levels falls above 4 mU/L but free T4 levels falls within the reference range; the state of the axis is said to be clinically subclinical hypothyroidism. Suppose TSH levels are above 4 mU/L and free T4 levels below 7 pg/mL; the state of the axis is said to be clinically hypothyroidism (an underactive thyroid gland). Keeping free T4 levels within the normal reference range is very important . Lower levels of free T4 can cause hypothyroidism that results in several health problems including obesity, joint pain, infertility, slow metabolism, puffy face, constipation, stiffness, dry skin, depression, fatigue and higher heart rate .\n
Hyperthyroidism is another clinical state of the axis, in which free T4 levels stay above 18 pg/mL while the levels of TSH stay below the reference range 0.4 mU/L . This state of the axis occurs when the thyroid gland over produces and secretes the thyroid hormones either in response to chronic TSH stimulation or due to Graves’ disease. Transient hyperthyroidism is a phenomenon that refers to the leakage of stored thyroid hormones from the gland, typically called Hashitoxicosis or thyroid burst . It happens due to the Hashimoto’s autoimmune thyroiditis. The bursting of thyroid gland can happen at any clinical stage of the patients.\n
Hashimoto’s autoimmune thyroiditis is one of the immune disorders hosted by the thyroid gland [6, 7]. Under this disease process, the thyroid follicular tissue undergoes the slow destruction by the immune system and thereby the thyroid struggles to produce enough hormones for the body’s requirement. Currently, the incidence rate of Hashimoto’s autoimmune thyroiditis is estimated to be 300–500 cases per 100,000 individuals per year. In this disease, the interaction of the immune system with the thyroid is highly complex involving the cellular and humoral mechanisms. The involvement of humoral mechanism can be verified through a simple blood test for anti-thyroid antibody biomarkers. More specifically, the elevated titers of anti-thyroid antibodies, thyroid peroxidase antibodies (TPOAb) and/or thyroglobulin antibodies (TGAb) may indicate the malfunction of the immune system . For simplicity in patient-specific modeling, we choose TPOAb as a biomarker representing the presence of the Hashimoto’s autoimmune thyroiditis.\n
In 1912, Hakaru Hashimoto published the description of four cases of diffused goiter which were showing the signs of infiltration of immune cells and antibodies in the thyroid gland . All these goiters appeared differently from the colloid goiters that occur due to insufficient iodine in-take from the diet. In his time, the idea of autoimmune disease had not been established in the medical literature. He isolated these cases and asked the future clinicians to explore the mechanism of this new form of goiter [10, 11]. In 1936, the description of goiter due to lymphocytic thyroiditis was rediscovered in the United States and was labeled as Hashimoto’s thyroiditis. It has been characterized as an organ-specific autoimmune disease and swelling might be the result of the chronic stimulation of the thyroid gland by the serum TSH [12, 13]. Now, scientists have described autoimmune thyroiditis with atrophy, a variance from the original description of Hakaru Hashimoto.\n
The pituitary gland has no clue about the undergoing immune process in the thyroid gland and does not know how to adapt to the changing environment and the abnormal function of the thyroid [14, 15]. It simply tells the thyroid to keep up with the release of thyroid hormones for its signal. By responding through bulging its size, the poor thyroid gland needs to accomplish the task of producing enough levels of hormones if it could. What initiates the immune process against the thyroid gland? The answer is still largely unknown to the researchers, but they all speculate the immune process might be initiated through a complex combination of the genetic and environmental pollution . Herein, we care about the functional size of thyroid gland as opposed to the exact size of the thyroid gland, which is treated as a hidden compartment . The exact size of the thyroid gland includes the size of the parathyroid glands, the isthmus, capillaries, blood flow and so on whereas the functional size includes the follicular cells that contains thyroglobulin molecules in which hormones are stored. Basically, the functional size of the gland counts the part of the gland that is able to make and secrete hormones; obviously this information is not available in the clinical setting.\n
A measurement of TPOAb titers is required from the blood test to diagnose autoimmune thyroiditis and afterward, TPOAb titers is not usually measured as its levels not helpful in the treatment. As mentioned above, the axis function is determined from one-time measurement of TSH and free T4, respectively. However, the functional size of the thyroid gland cannot be measured through the laboratory experiments or with any medical tools, which does not exist in practice at least so far. Even if it exists, then it would be very expensive in terms of cost and time. However, the mathematical model can help by telling us the functional size from the clinical measurements of specific patient and from the physiological point of view. Moreover, the functional size of the gland is different for everyone, such as boy’s thyroid size versus men’s thyroid size or girl’s thyroid size versus women’s thyroid size and so on. Using the model and the clinical measurements of TSH and free T4, we will determine the initial functional size of the gland for a given patient. Goiter due to this disease occurs more frequently in women than men . Using the model, we describe the diffused goiter as the functional size increases to keep up with the normal production of hormones.\n
In reality, it is hard to perform experiments on patients in the laboratory setting and collect data to understand the underlying dynamics of the pituitary-thyroid axis in autoimmune problem. Hypothetically speaking, if the experiment is possible, then one needs to consider patients’ time, cost, safety, and variability. Two patients with similar characteristics from societies might react to an experiment in an unusual way [18, 19]. The inter and intra variability among patients causes a major problem for treatment and challenges the obtained information from labs. Using the model, we can carry out an experiment for a specific patient and the hidden dynamics of the axis can be investigated at the microscopic level for shorter or longer time-period. In general, the model is implemented through a computer program and the parameters are tuned to a specific patient and then the experiments can be performed with the test hypothesis . Also, we can identify sensitive and insensitive parameters that are responsible for the physiology of the axis under autoimmune thyroiditis, that causes the hypothyroidism and goiter. This is an effective and modern way to explore the complex interaction of the immune system to the thyroid and its consequence on the negative feedback loop [21, 22].\n
The hormone TSH has a half-life of 1 hour in serum on the average scale. Similarly, the hormone free T4 has a half-life of 7 days in serum on average, the thyroid per-oxidase antibodies (TPOAb) has a half-life of 24 hours in serum on average and the half-life of functional size of the thyroid gland is not known and probably varies among individuals. Several time scales have been involved in the disease process. Based on the fundamental principle of rate laws, the model will be developed here to predict the history of patient-specific dynamics of the axis and thyroid size in autoimmune thyroiditis [23, 24, 25]. The model can unfold the consequences of the presence of antibodies titers in the blood. Basically, it can replace humans but mimics the system that provides convenience and flexibility for scientists to run experiments on the computer. Information obtained through the dynamics might be useful in treating patients and improving the accuracy of the levothyroxine treatment . For instance, the levothyroxine drug can be targeted and administered to the specific patient so that the levels of TSH and free T4 maintained within the normal reference range, which in turn can avoid other health problems.\n
The remaining of this chapter is divided into three main sections. Section 2 introduces the development and construction of a coupled model describing the interactions of two subsystems (the axis and humoral immune system) and reduces the coupled model based on the physiological assumptions. Section 3 analyzes the solutions of the reduced model qualitatively and numerically for various parameter values. Section 4 describes the case studies of three patients, validity of the model and provides the conclusion of this work.\n
The pituitary gland is diseased free, so the feed forward loop is intact.
Total TSH receptors concentration does not change during Hashimoto’s autoimmune thyroiditis.
Serum TSH stimulates the growth of functional thyroid and the production and secretion of thyroid hormones.
The humoral immune system uses serum TPOAb to attack the thyroid and those titers can be used as a biomarker for the level of the anti-thyroid immune activity.
The patient does not demonstrate central or peripheral resistance to thyroid hormone.
Serum free T4 have much slower dynamics in the disease process compared to serum TSH, serum TPOAb and the functional size of the gland.
We denote \n\n, \n\n and \n\n to represent the concentrations of serum TSH (mU/L), serum free T4 (pg/mL) and serum TPOAb (U/mL), respectively. We let \n\n represent the functional size of the gland. According to the fundamental principle of rate laws, the rate of change of concentration of TSH over time \n\n is the secretion rate from the pituitary gland minus the elimination rate of TSH via kidney through the unspecified mechanism. The rate of change of concentration of free T4 over time \n\n is the secretion rate of free T4 from the thyroid minus the elimination rate of free T4 via kidney through unspecified mechanism. Similarly, the rate of change of concentration of TPOAb over time \n\n is the production rate due to the abnormal interaction thyroid gland and immune system minus the elimination rate of TPOAb via kidney through unspecified mechanism. Finally, the rate of change of the functional size of thyroid gland over time \n\n is the growth rate minus the destruction rate of the gland.\n
Next, we write a coupled model that represents the interaction of two subsystems: the negative feedback loop and the humoral immune system.\n
where \n\n and \n\n. The normal average value for each state variable is given in \nTable 1\n. The model has 11 positive parameters unchanged over time, in which each parameter has a unique physiological meaning and their corresponding values listed in \nTable 2\n. The model (1)–(4) captured the dynamics of two subsystems coupled through the functional size of the thyroid gland. Using assumption (8) that serum free T4 has much slower dynamics compared to serum TSH, serum TPOAb and the functional size of the thyroid gland, we can take \n\n and obtain a reduced model consisting of three differential equations and one algebraic equation:\n
|Name\n||Normal value\n||Normal range\n||Source\n||Unit\n|
|\n\n\n\n\n||\n\n\n\n\n||\n\n\n\n\n||Literature \n||mU/L day\n|
|\n\n\n\n\n||\n\n\n\n\n||\n\n\n\n\n||Simulation\n||Pg/mL L day\n|
Using (Eq. (6)), the functional size can be determined for patients if a result of TSH (\n\n) and free T4 (\n\n) value is known from the blood test, which in turn can be used as an initial condition on (Eq. (3)) for the model simulation. The graph of (Eq. (6)) is shown on the \nFigure 2\n, which illustrates the functional size of thyroid gland in terms of varying TSH and free T4 values. Using the reduced model, we can describe the clinical progression of overt hypothyroidism and thyroid size such as goiter and atrophy.\n\n
which leads to diseased-free and diseased steady state solutions besides the initial condition. The first steady state (diseased-free) denoted as \n\n where \n\n and the value of \n\n is given by the cubic equation:\n
By Descarte’s rule of signs, Eq. (7) has one positive real solution, so the reduced model has diseased-free state in the positive octant for the system parameters. In fact, this steady state solution lives on the surface of the function \n\n given by Eq. (6). The operation of the pituitary-thyroid axis is healthy and fully functional in the presence of diseased-free steady state solution.\n
Next, the second steady state is the diseased state solution denoted as \n\n where\n
and the value \n\n is given by the quadratic equation:\n
By Descarte’s rule of signs, Eq. (8) has one positive real solution when \n\n so the reduced model has diseased state in the positive octant for the system parameters.\n
We let \n\n be the unique value of N, where\n
the system undergoes a bifurcation. We call \n\n is a bifurcation value of the system (Eqs. (3)–(5)) and \n\n is the bifurcation parameter. We can take the initial state of the system at where the diseased-free and diseased steady states merge together. Using this definition, we can also calculate the parameter \n\n (say \n\n) provided the values of \n\n and \n\n is known.\n
We first define the level surface (S) through a function of three variables from Eq. (6) as\n
The normal surface vector \n\n can be calculated from the gradient of this function\n
In particular, the normal vector of the surface \n\n at initial state \n\n is a vector perpendicular to the tangent plane of \n\n at \n\n is given by\n
The equation of the tangent plane at point \n\n is given by\n
Next, we take the implicit differentiation of Eq. (9) with respect to time \n\n and substituting the initial point \n\n, which results in the following equation for \n\n:\n
Notice that \n\n is independent of the parameters \n\n and \n\n. So, we can use this value of \n\n as \n\n. This threshold value can be uniquely determined for every patient based on their initial blood test results and if we know other parameters \n\n, \n\n and \n\n. We have estimated the values of all other parameters from the literature or through calculation (see \nTable 2\n).\n
As the bifurcation parameter \n\n varies from the threshold value \n\n, the number of the steady states changes, and its behavior near the steady states (linear stability) changes as well. Moreover, as the bifurcation parameter varies, the appearance and disappearance of steady states can be seen on the tangent plane at the initial state of the function \n\n. For this system (Eqs. (3)–(5)), there are three cases, namely when \n\n, \n\n and \n\n. \nFigure 3\n illustrates a nice overview of how the solutions become stable and unstable as parameter \n\n changes from 0 to 300.\n\n
When \n\n, the system (Eqs. (3)–(5)) has two steady states, which undergoes a saddle-node bifurcation—disappearance of the saddle diseased state. One can imagine if we start the system from the initial steady state, which splits into the diseased-free and diseased state. The disease state is asymptotically stable which means all its eigenvalues are purely negative real. However, the diseased-free state is saddle because two of its eigenvalues are positive real and one of its eigenvalues is a negative real. See  for proof of the stability of the diseased-free and diseased steady state.\n
As the value of \n\n decreases from the threshold value and when it is close to zero, there exists a homoclinic orbit associated at the initial state which forms the boundary for the asymptomatically stable interior diseased state. The orbit captures the hidden dynamics of the diffused goiter and hashitoxicosis (see \nFigures 4\n and \n5\n).\n\n\n
When \n\n, the system (Eqs. (3)–(5)) has one steady state (diseased-free), which is asymptotically stable. Since there is only one steady state in the system, the operation of the pituitary-thyroid axis in the presence of anti-thyroid peroxidase antibodies can be tested. The numerical simulation of the model shows TPOAb slowly approaches zero as time increases over months (see \nFigure 6\n).\n\n
When \n\n, the system (Eqs. (3)–(5)) has one steady state (diseased-free), which is asymptotically stable. As the value of \n\n increases beyond the threshold value, the diseased free state is moving on the function \n\n describing the progression of thyroid size and the function of the axis. One can imagine if we start the system from the initial state of normal thyroid size, then the system may approach goiter or atrophy depending upon the parameter value of \n\n while keeping other parameters fixed. On the other hand, varying the parameter \n\n, the clinical mechanisms euthyroidism \n\n subclinical hypothyroidism \n\n overt hypothyroidism can be explained while keeping \n\n fixed . See \nFigure 7\n that shows the clinical progression of subclinical hypothyroidism and thyroid size.\n\n
Suppose that an individual is diagnosed with Hashimoto’s autoimmune thyroiditis; we can use the reduced model to explain the physical and clinical symptoms occurring for this individual in the course of this disease. Having this disease means everyone has a unique behavior and knowing that may be very helpful in managing the course. In fact, certain parameter values in the reduced model are responsible for the uniqueness in patient behavior. More specifically, we have identified through the stability analysis that the parameter \n\n can explain the size of the thyroid (goiter, normal or atrophy) whereas \n\n can explain the clinical progression of the disease [21, 22]. Putting these two parameters together, one can test the hypothesis on the patients’ symptoms. See \nFigure 8\n that shows the parameter curve in terms of \n\n and \n\n that is responsible for different thyroid size and clinical conditions.\n\n
Suppose a patient’s information is given; then a threshold \n\n can be found for the patient using Eq. (10) and the definition of the bifurcation parameter. If another information is available in the future, then another threshold can be obtained. To be exact, the threshold is an ordered pair found on the curve given by the bifurcation definition of parameter (see \nFigure 8\n). To simulate physical and clinical symptoms, we take two test ordered pairs from the parameter curve below and above the threshold \n\n) while other parameters came from \nTable 2\n. With the first test point for \n\n, we simulated the model for a day and found the symptoms of the disease are hashitoxicosis and goiter (see \nFigure 9\n). Similarly, by keeping the second test point as (5.8, 397.9), the one-day simulation showed mild atrophy and overt hypothyroidism (see \nFigure 10\n). Finally, by keeping the third test point away from the curve as (4.5, 150), the one-day simulation showed subclinical hypothyroidism and normal size of the thyroid (\nFigure 11\n).\n\n\n\n
Using patients’ information from the peer-reviewed published article, we will predict patients’ natural history of the disease. More precisely, the reduced model can describe thyroid size and clinical progression from euthyroidism to subclinical or overt hypothyroidism for each patient given below. These patients’ information was provided by Salvatore Benvenga originally and were already published in the article . The data consists of TSH, free T4 and TPOAb information whose normal reference ranges are (0.4–2.5) mU/L, (7–18) pg./mL and (0–200) U/mL, respectively. Using Eq. (6) and data, we have computed the functional size of the thyroid gland. Using Eq. (10), data and the bifurcation definition, we have computed \n\n and \n\n for all patients given below (see \nTables 3\n–\n5\n).\n
|Time (months)\n||TSH (mU/L)\n||FT4 (pg/mL)\n||TPOAb (U/mL)\n||Thyroid size (mL)\n||\n\n\n\n\n||\n\n\n\n\n|
|Time (months)\n||TSH (mU/L)\n||FT4 (pg/mL)\n||TPOAb (U/mL)\n||Thyroid size (mL)\n||\n\n\n\n\n||\n\n\n\n\n|
|Time (months)\n||TSH (mU/L)\n||FT4 (pg/mL)\n||TPOAb (U/mL)\n||Thyroid size (mL)\n||\n\n\n\n\n||\n\n\n\n\n|
Assuming all patients had diseased-free state (1, 13, 0.015, 0) at some point, the data provides us the diseased state from which one can back trace the physical and clinical conditions caused by the disease and experienced by the patients. We have taken three untreated patients with a known biomarker for the presence of the autoimmune thyroiditis. Patients 1 and 2 visited the clinic three times whereas patient 3 visited the clinic six times due to mild goiter or other symptoms. We do not have any information about their thyroid sizes in the data.\n
In \nTable 3\n, the clinical measurements, the functional thyroid size and parameters are listed for patient 1 which form the diseased states. The model simulation started from each diseased state for 10 years to back trace the natural symptoms of patient 1 due to autoimmune thyroiditis. The course of disease revealed thyroid size changes from mild goiter to normal while the function of axis remained normal (see \nFigure 12\n). In \nTable 4\n, the clinical measurements, the functional thyroid size and parameters are listed for patient 2. The model simulation is done for 1 year to back trace the symptoms from the course of the disease. The first disease course revealed a development of subclinical hypothyroidism with normal functional thyroid size (see \nFigure 13\n). Similarly, patient 3 diseased states are found in \nTable 5\n. The model simulation is done only for 2 days to capture the early course of the disease. Surprising results are seen in the early course of the disease (see \nFigure 14\n). This patient had developed both clinical and physical symptoms at different points in time. In particular, the first course shows mild goiter and mild hashitoxicosis whereas the sixth course shows the overt hypothyroidism and mild atrophy.\n\n\n\n
The human body is made up of so many subsystems such as the pituitary-thyroid axis and the immune system. These subsystems do not disrupt the function of each other in healthy people. A function of the pituitary-thyroid axis is to secrete the appropriate levels of thyroid hormones and take care of body’s metabolism whereas the function of the immune system is to protect the body’s organs such as the thyroid gland and remove foreign substances that try to invade the organs. Hashimoto discovered the abnormal interaction of the immune system to the thyroid gland, which resulted in the disruption of the physiology of the axis. Later this diseased condition has been named as Hashimoto’s autoimmune thyroiditis. Hashimoto noticed the destruction of thyroid follicular cells through physical and clinical symptoms in four middle aged female patients.\n
A keystone of the functional thyroid gland is the follicular cells, which die due to the aggressive and destructive attack by the immune system. Hashimoto patients discovered the disease sometimes by themselves due to discomfort in the neck (small goiter) or accidentally during annual checkup by the family physicians. Consequences of Hashimoto disease can be classified into physical and clinical symptoms at various stages of the disease. Both types of symptoms occur in a sequential manner from one extreme to another. More precisely, the physical symptom runs from goiter \n\n atrophy whereas the clinical symptom runs from hashitoxicosis \n\n overt hypothyroidism. In order to describe these symptoms mathematically, we chose four state variables namely the size of functional thyroid gland to describe the physical symptom, TSH and free T4 to describe the clinical symptoms and thyroid peroxidase antibodies as a representative for the presence of Hashimoto disease. Modeling the disruption of the axis through a coupled model is the key to unlock the hidden dynamics experienced by the patients. The hidden dynamics can be seen in patients through physical and clinical symptoms.\n
The literature description of Hashimoto disease begins with a gradual swelling of the thyroid gland and development of mild clinical condition, euthyroidism or subclinical hypothyroidism and subsequent gradual progression of overt hypothyroidism. Small goiter and hashitoxicosis are the very early stages of the disease and typically go untreated and hidden in the view of patients and physicians. Overt hypothyroidism is the irreversible end clinical state (where levothyroxine treatment is needed) whereas atrophy is irreversible end physical state of the disease. Basically, the mechanism involved in the progression of the disease is unique and sequential. For instance, some patients may have the disease course of goiter and euthyroidism, some may have goiter and the clinical progression euthyroidism \n\n subclinical hypothyroidism and some patients may have gradual progression from normal thyroid size to atrophy and euthyroidism to overt hypothyroidism.\n
Herein, we have developed and used patient-specific model to describe all possible mechanisms involved in the autoimmune thyroiditis. This can be achieved using two parameters \n\n and \n\n. To be precise, the parameter \n\n and \n\n can trace the thyroid size and clinical progression. We validated the model with three patients’ data by describing their natural history of the disease.\n
I would like to dedicate this chapter to Dr. Stephen J. Merrill for his guidance on mathematical modeling of Hashimoto’s autoimmune thyroiditis.\n
We declare there is no conflict of interest on this chapter.
I would like to thank Dr. Salvatore Benvenga for providing Hashimoto’s patients data for the modeling work.\n