Mathematical-based morphological classification of skin eruptions corresponding to the pathophysiological state of chronic spontaneous urticaria

Background Chronic spontaneous urticaria (CSU) is one of the most intractable human-specific skin diseases. However, as no experimental animal model exists, the mechanism underlying disease pathogenesis in vivo remains unclear, making the establishment of a curative treatment challenging. Methods A novel approach combining mathematical modelling, in vitro experiments and clinical data analysis was used to infer the pathological state of CSU patients from geometric features of the skin eruptions. Results Based on our hierarchical mathematical modelling, the eruptions of CSU were classified into five categories, each with distinct histamine, basophils, mast cells and coagulation factors network signatures. The analysis of 105 real CSU patients with this classification by six individual dermatologists achieved 87.6% agreement. Furthermore, our network analysis revealed that the coagulation status likely determines boundary/area pattern of wheals, while the state of spontaneous histamine release from mast cells may contribute to the divergence of size and outline of the eruptions. Conclusions Our multi-faceted approach was accurate in defining pathophysiological states of disease based on geometric features offering the potential to improve the accuracy of CSU diagnosis and better management of the disease in the clinic.


Plain language summary
Chronic spontaneous urticaria (CSU) is a persistent skin disease that causes red itchy skin eruptions, called wheals, that are of various shapes.These wheals repeatedly appear and disappear daily for up to weeks or even decades, severely impacting the quality of life of those affected.The causes, consequences or disease processes are largely understudied.Here, we developed a novel approach using mathematical modelling to analyse the geometric measurements of patients' wheals alongside laboratory experiments on human skin samples with eruptions to improve CSU diagnosis and clinical management.We find that geometric measurements of these skin eruptions can be classified into five categories, which could facilitate accurate diagnosis and consequently better clinical management of CSU and may be more widely applied to other skin diseases.
T he skin is the largest organ in the human body and plays an important role in maintaining homeostasis and protecting the body from the external environment.Antigens, chemicals, and ultraviolet rays, and endogenous substances can cause cutaneous inflammatory reactions and tumorigenesis.A variety of skin diseases have an impact on society due to seriously impairing the quality of life of patients, even if they are not life-threatening.Urticaria is one of the most common human skin disorders, affecting at least one in five people during their lifetime 1,2 .It is characterized by the appearance of skin eruptions, also known as wheals, that are of various shapes ranging from a few millimeters to several centimeters in size and is typically accompanied by itching.Patients with urticaria suffer from these symptoms potentially for years to even decades.Among the different subtypes of urticaria, chronic spontaneous urticaria (CSU) is characterized by the spontaneous appearance and disappearance of wheals 3 .
The pathophysiological states associated with CSU include autoimmune responses, cellular infiltrates, coagulation, and the complement system 4 .However, how these responses contribute to wheas remains unknown, which is largely due to the limited understanding and scarcity of clinical data on the in vivo pathological structure of urticaria.Furthermore, due to the human-specific nature of the disease, there is also a lack of an appropriate animal model system that can recapitulate the entire disease pathology 5 .Consequently, elucidating the underlying pathophysiological dynamics in vivo remains challenging, and developing patient-specific effective treatments for urticaria remains difficult.
To overcome the limitations of the lack of human-specific skin disease models, we developed a novel approach linking the geometric features of skin eruption with pathophysiological structures in vivo.As all skin diseases are characterized by "visible outputs", skin eruptions can be considered a reflection of in vivo pathophysiological dynamics.However, eruption geometry has neither been fully understood nor explored as a measure of the differences in pathological states between patients.
Thus, in this study, we focus on the morphological features of skin eruptions and develop a novel approach linking the eruption geometry of patients to the in vivo pathological dynamics of CSU using hierarchical mathematical modeling, in vitro experiments, and clinical data analysis of the eruption patterns of urticaria.We show that the major in vivo pathophysiological state of CSU can be typically defined as one of five skin eruption pattern types, thereby predicting degrees of contribution by each pathogenetic factor in urticaria.This study suggests that mathematical modeling analysis of morphology appearing on the body may help not only with diagnosis and classification but also reveal molecular mechanisms and targets for new treatments of a variety of human diseases that develop geometric lesions.

Methods
Collection of images of wheals and their morphological classification.Images of skin eruptions (wheals) were collected from 105 patients with CSU who visited the Department of Dermatology at Hiroshima University Hospital between January 1, 2000 and April 15, 2022.The institutional ethics review board approved the study protocol (Ethical Committee for Epidemiology of Hiroshima University, Hiroshima, Japan, approval number; E2020-2388), which involved analyzing images stored in the secured hospital database.Informed consent was obtained in the form of opt-out on the website.Six dermatologists, specializing in allergic skin diseases, classified the morphology of wheals on each CSU patient into one of the following seven major categories: annular, broken-annular, geographic, circular, dot, uniform (shapeless and spread throughout the body), and none of the above (NA).
The uniform and NA cases were considered unclassifiable based on the EGe criteria; therefore, six categories were used in the final analysis.Sample number 105 exceeded the minimal number required for an error with the true classification result within 8% and the significance level (α) was 5% 6 .
Statistical analysis.To classify the eruption patterns of patients, we used the vote-counting method, performed by six dermatologists.The reliability SB was calculated as follows: SB is the probability of pðA ^BÞ=pðAÞ where p(A) and p(B) are the probabilities that eruption patterns are classified in one of two classes or one of five types by doctors, respectively (Table 1).The reliability SS was calculated as follows: SS is the probability that both the class and type of eruption pattern were chosen correctly and that the SS-based analysis would provide more reliable statistical results.To confirm the statistical verification of the classifiable case data, we calculated the P-value as follows: where k 0 is the sample number of classifiable cases, n is the total sample number, and p 0 is the probability of the null/alternative hypothesis that the classifiable rate in the population is below/ above p 0 %.
Function representing the skin eruption state.To reflect the eruption state of the skin ð e ΩÞ based on the molecular dynamics of the dermis ðΩ D Þ, we define the eruption state function, S w : Ω D ! e Ωð&R 2 Þ, as follows: This function connected the dynamics of the coagulation factors ( A ½ ðx;tÞ ½Cðx; tÞ) leaked from blood vessels or histamine ( A ½ ðx;tÞ ½H M ðx; tÞ) released from mast cells in the dermis to the state of wheal formation seen on the skin ( e Ω). β w is the positive value and A ½ r is the threshold value at which the eruption appears on the skin.This function allows for two states of eruption states: eruption present (1) or absent (0) 7 .Because the patterning dynamics of histamine concentration of mast cells and coagulation factors were identical (Supplementary Figs.S2-S6), we plotted the representative results as A ½ ðx;tÞ ½Cðx;tÞ in this study.
For the impedance analysis, HUVECs were harvested using trypsin and seeded onto E-Plates (Roche Applied Science, Upper Bavaria, Germany) at a density of 50,000 cells/well.On day 2, cells were treated with LPS (100 ng/mL) and histamine (10 μM) only for the TF stimulation assay.On day 3, the E-plate was set on an iCElligence microelectronic biosensor system (Roche Applied Science) and impedance termed "cell index (CI)" was measured every 10 s.For the TF stimulation assay, cells were stimulated with indicated TF concentrations in DMEM/F12 with 2% BSA and 1% factor VIII-deficient human plasma without any supplements and FCS.For the histamine stimulation assay, cells were stimulated with indicated concentrations of histamine in DMEM/F12 supplemented with 10% FCS.
Statistics and reproducibility.Inhibition ratios of histamine release from human mast cells and basophils are shown as means with standard errors of three independent experiments.Gap formation of endothelial cells analyzed by impedance sensors is shown as a representative of three independent experiments.Statistical analysis was performed using Microsoft Excel TM (Version 16.78) in this study.
Reporting summary.Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.

Results
CSU mathematical modeling based on pathophysiological structures.To design our CSU mathematical model, we first defined all the necessary pathophysiological molecules and their interactions that have been associated with CSU, based on in vitro experimental studies (Fig. 1).Firstly, we took into consideration that external stimuli or internal auto-allergens and auto-antibodies stimulate both basophils and mast cells to release mediators by crosslinking the high-affinity IgE receptors (FcεRI) on their surfaces [8][9][10] [Fig.1a (0)].Next, this process subsequently enhances TF expression on vascular endothelial cells [Fig.1a (1)], which induces the extrinsic coagulation cascade [11][12][13] [Fig.1a (2 a )].Consequently, the coagulation factors increase the degranulation of basophils [Fig.1a (2 b )] or skin mast cells via gap formation by vascular endothelial cells [Fig.1a (2 a ), (3), (4)] leading to edema [Fig.1a (5), (6)] 14 .Furthermore, adenosine, which is simultaneously released following degranulation of both human peripheral basophils and skin mast cells 15 , suppresses histamine release and TF expression 16 , we thus also took this into consideration (Fig. 1a, green inhibition symbols).Based on these facts, we applied the pathophysiological structure into a network structure (Fig. 1b), focusing on the following main factors that are expected to have a key role in the in vivo network: (i) concentration of histamine released from basophils (½H B ), (ii) TF expression on vascular endothelial cells ( TF ½ ), (iii) coagulation factors leaked from blood vessels ( C ½ ), and (iv) histamine released from mast cells (½H M ).We developed the hierarchical mathematical model integrating both the intravascular and extravascular dynamics using in vitro Annular pattern Annular eruptions surrounded by a single arc with a sharp margin.They may be fused with each other.

Broken-annular pattern
Petal-like annular eruptions composed of multiple arcs.There is no continuity or intersection of individual arcs expanding centrifugally.

Area pattern Inflated eruptions with/ without flatted uniform erythema
Geographic pattern Circular shaped wheals with a tendency to fuse with each other, forming a geographic appearance.
Circular pattern Isolated circular wheals of ≥1 cm diameter with a low tendency to fuse together.

Dot pattern
Small wheals or erythema of <1 cm in diameter, with a low tendency to fuse together.
Morphologically, eruptions in chronic spontaneous urticaria (CSU) are divided into two classes and five subtypes.
experimental data.Wheals in CSU are caused by dilatation and increased permeability of the postcapillary venules in the subpapillary layer of the dermis.The subpapillary layer of the dermis is the particularly shallow layer of the dermis with an average of 2 mm thickness, which depends on location.
Compared with the scale that reproduces the horizontal spread of the wheals in the mathematical model (35 cm × 35 cm), the depth scale at which the vascular reactions involved in wheal formation is negligibly small.Thus, we mathematically assumed that the dermis of the skin can be considered as a thin twodimensional (2D) sheet with numerous capillaries spread over it.
We constructed a 2D spatial model with the intravascular space distributed homogeneously throughout the dermis.Based on the network (Fig. 1b), we constructed the following model (see Supplementary Methods 1-5 for details on model construction and parameter estimation).
where Ω D ð& R 2 Þ is the dermal region, Ω E is the vascular endothelial tissue, and Ω B ð& Ω E Þ is a randomly chosen local region of vascular endothelial tissue where the basophils adhere to the endothelial cells.Furthermore, δ B ; δ T ; and δ M are the basal production rates and μ B ; μ T , and μ M are the basal decay rates.γ B and γ M are histamine release rates of basophils and mast cells, respectively (network ( 2) and ( 4) in Fig. 1a, b).The terms χ B x; t ð Þ and χ M x; t ð Þ are the Heaviside step functions, assumed to be 1 when each basophil and mast cell contained histamine and 0 when they each released their total histamine content 7,17 .We quantitatively estimated the main network functions, g B inhibition ; g T inhibition ; g M inhibition ; f T activation , and f leakage x; t ð Þ based on the in vitro experiment data.g B inhibition ; g T inhibition ; and g M inhibition represent the inhibition rates of adenosine which is produced by adenosine triphosphate (ATP) released simultaneously with and in proportion to histamine from basophils and mast cells 15 (green inhibition symbols in the networks shown in Fig. 1a, b).Adenosine is known to suppresses the histamine release of human peripheral basophils and human skin mast cells (hsMCs) in response to anti-IgE in a concentration-dependent manner 16 .Based on our previous data 16 , we directly analyzed the inhibitory rates of histamine released from hsMCs and basophils via adenosine (Fig. 2a, b, respectively).Using these data, we estimated the specific forms of g B inhibition and g M inhibition as shown in Fig. 2c: where α B and α M are the maximal inhibition rates, α B0 and α M0 are the constants that determine the degree of increase in the curves.
Next, using in vitro experiment data of the tissue factor (TF) suppression by adenosine 18 , we also analyzed the inhibitory rates of TF mRNA expression based on adenosine concentrations (Fig. 2d) and estimated the inhibitory rate function of TF follows (Fig. 2e): where we estimated the parameter values, the maximal inhibition rate (α T ) and a constant ðα T0 Þ, by the least squares method.f T activation represents a function of histamine-induced TF expression (networks (1) in Fig. 1a, b).TF in human umbilical vein endothelial cells can be induced by histamine and synergistically increased by co-stimulation with histamine and lipopolysaccharide (LPS) 18 .To reflect the potential synergistical output of TF activation in our model, we directly analyzed previous experimental data on the relative expression level of TF mRNA induced by co-stimulation with histamine and LPS (Fig. 2f).We then converted the result to data based on the histamine-only concentration and estimated the TF activation function as follows (Fig. 2g): where γ T is the maximal activation rate and γ T0 is a constant that determines the degree of increase in the curve.
The TF expressed on endothelial cells induces the extrinsic coagulation cascade 11,13 , and activated coagulation factors produced by the extrinsic coagulation pathway induce intracellular gap formation between endothelial cells via PAR-1 19 .Thus, the coagulation factors leaked from intravascular, namely, f leakage x; t ð Þ can be determined by the dynamics of intercellular gap formation.
In order to determine the detailed form of f leakage x; t ð Þ, we performed in vitro experiments using human umbilical vein endothelial cells (HUVECs).We looked for the changes in the adhesion area and morphology of the cells cultured on the electrodes with varying TF and histamine concentrations, respectively.We first measured the time-related cell index (CI) with varying TF by using an impedance analysis (Fig. 2h).We found that the CI was not significantly changed following treatment with TF at <1 ng/mL, whereas it was dramatically decreased with a TF concentration >3 ng/mL.Notably, the change in CI was classified at two levels and its dynamics was switched by a threshold concentration of TF between 1 and 3 ng/mL.This result suggests that TF plays a role in inducing the leakage of coagulation factors via switch-like gap formation.
Based on this experimental result, we inferred the TF-dependent switch-like function, χ sw TF ½ ð Þ; by calculating the relative CI, ð1 À CIÞ=ð1 À CI 100ng=ml Þ, which gives the relative variation in the size of the CI for each time point.We calculated the average relative CI at [0.2, 0.6] (hour) for each TF concentration and found the switch-like function as shown in Fig. 2I, j.Furthermore, we found that the rate of gap area expansion is determined in a concentration-independent manner for both TF and histamine (Supplementary Method 5, Fig. S1 in mathematical and experimental details).Taken together, we obtained where γ C is a leaked rate of coagulation factors (Fig. 2j).
At this stage, we had completed the CSU mathematical model ( 5)-( 8) based on experimental analysis of the key pathophysiological networks.The eruptions of CSU visibly appeared on the skin where the leaked coagulation factors from blood vessels or histamine concentration were high in the dermis.Therefore, we plotted the eruption state of the skin using the eruption state function (defined in the Methods).
Recapitulation of eruption patterns in urticaria using a mathematical model.We next validated whether our mathematical model could successfully recapitulate the CSU eruption patterns.Depending on the parameter values, we identified several eruption patterns observed in the patients with urticaria (Fig. 3a-e, Figs.S2-S6).Based on geometric characteristics, we categorized the different eruption patterns identified in our mathematical model as follows: annular, broken-annular, geographic, circular, and dot patterns.
The eruption image data of patients in clinical medicine represent the data of a certain time point of eruption development dynamics.Thus, we compared the in silico eruption patterns at some fixed time point with the image data collected for CSU patients who visited the hospital (Fig. 3f).We found that patterns among the patient data were highly similar to all pattern types observed in the in silico experiments, hence validating our CSU mathematical model.

Classification of eruption geometry in urticaria patients
according to eruption geometry criteria.Interestingly, during the in silico validation of our mathematical model, we observed that eruption patterns were typically identified as one of five types when parameter values were varied (Fig. 3g, h).Therefore, we next sought to determine whether all eruption patterns appearing in CSU patients could be classified into specific types and then linked to in vivo pathologies.To this end, we first developed classification criteria defining the clinical features of the eruption (wheal) geometry found in the in silico patterns (Table 1).We designated the resultant table of clinical features Criteria for Classification of Eruption Geometry (EGe criteria).We classified the patient eruption patterns into two large classes, boundary and area patterns, which we then further classified into five specific patterns, namely annular, broken-annular, geographic, circular, and dot patterns, as observed in the in silico experiments.
To confirm whether the eruption patterns of actual CSU patients could be classified into one of the eruption types in the EGe criteria, we collected the eruption image data of 105 patients who visited the Hiroshima University Hospital with CSU-related skin eruptions (see Online Method).The data were then analyzed by six dermatologists, based on the EGe criteria.We evaluated the reliability SB (≡[number of doctors who selected the sample as one of the types in the EGe criteria]/[number of dermatologists ( 6)]), of the statistical result using the vote-counting method.
The results showed that 87.6% of the samples could be classified as one of the eruption types with 100% SB reliability (Fig. 4a), and 0% of the samples were unclassifiable, with a reliability SB > 83% (Fig. 4b).We also statistically confirmed that the classifiable rate of 87.6% is not a random chance by the Pvalue (Supplementary Method 6).This result proves that the eruption types in CSU patients could be classified into one of the patterns of the EGe criteria developed based on our mathematical model.This suggests that the eruption dynamics of CSU patients can be determined using our CSU mathematical model.
Next, we investigated which types of eruption patterns were commonly observed in CSU patients.To obtain results with higher reliability, we introduced a new measurement of reliability (SS) defined as SS = [SB of type selection] × [SB of class selection] (Fig. 4c, Online Methods).From the analyzed data with the reliability SS, we found that geographic and circular types were observed at high frequencies, whereas boundary patterns appeared fewer times than area patterns did.
Linking the key features of eruption geometry to related pathological networks via a mathematical model.We next extracted the minimal sets of key characteristics of eruption patterns from the features defined in the EGe criteria (Table 2): distinguishable (KF1)/indistinguishable (KF3) between the boundary and interior of the wheal, disconnectedness in the wheal boundary (KF2), fusion of the wheals (KF4), and no change in size of wheals (KF5).This showed that each of the eruption patterns could be identified using appropriate combinations of the five features.Specifically, the boundary pattern is identified by KF1, the annular/broken-annular pattern is identified by KF1 and KF2, the geographic/circular pattern is identified by KF3 and KF4, and the dot pattern is identified by KF3 and KF5.
We then looked at how each eruption pattern is connected to the network of pathological dynamics.To this purpose, we performed a global sensitivity analysis of the Extended Fourier Amplitude Sensitivity Test (eFAST) ( 20 ; explained in detail in Supplementary Method 7), which tells us which networks play an important role for a given sensitivity function, namely, each eruption pattern.Based on the five key features of clinically observed wheals in CSU (Table 2), we defined five sensitivity functions to measure the degree of alteration of the features of each eruption pattern in response to changes in parameter values representing the network intensity (Fig. 5, Supplementary Method 7 shows the details of sensitivity functions).
We first noted that most key features had a relatively high association with the networks related to TF and mast cells (Fig. 5f, yellow-colored regions), suggesting that the dynamics of TF and mast cells may be important players in determining the geometrical properties of CSU eruptions.We next compared the difference between boundary and area features.In the boundary feature (KF1), the decay rate of TF ðμ T Þ highly influenced the features (Fig. 5a, f).By contrast, we found that μ T had a negligible effect on the area feature (KF3) (Fig. 5c, f).In addition, the activation network of TF (γ T and δ T ) had a considerable effect on the area feature but not on the boundary feature, although an effect was evident (Fig. 5f).These results suggest that the feature of boundary/area pattern is likely to be highly influenced by TF degradation/activation.
Next, to determine the network that played a critical role in the divergence of the annular and broken annual patterns, we compared the sensitive analysis result between the boundary (KF1) and boundary disconnectedness (KF2) features (Fig. 5a, b).We found that the effect of the basal activation network (δ M ) of histamine release from mast cells in the disconnectedness feature (KF2) was not as high as that of the boundary feature (KF1), indicating that the state of spontaneous histamine release from mast cells may be the critical factor to make the divergence of the annular and broken annual patterns.
We then compared the area fusion (KF4) and punctate (KF5) features (Fig. 5d, e).The first notable differences between these two features were the activating effect of spontaneous histamine release from mast cellsðδ M Þ and the inhibitory effect of histamine release from mast cellsðα M Þ.Both the networks (δ M and α M ) showed a high effect on the area fusion feature, whereas they exhibited relatively negligible effects on the punctate feature.This result suggests that the area types of eruption patterns may become geographic/circular type or dot type depending on whether the histamine dynamics from mast cells are strongly involved or not.We summarized the classification of eruption types and their relations with TF and histamine dynamics of mast cells based on the sensitive analysis in Fig. 5g.

Discussion
In this study, we established a novel approach to overcome the limitations of clinical data by developing a three-step research strategy: generating a mathematical model that recapitulates eruption patterns from the pathological structures (Fig. 5h, Step 1), establishing EGe criteria for the clinical classification of eruptions in patients based on in silico eruption patterns (Fig. 5h, Step 2), and linking the features of eruption geometry to critical networks of pathological structures (Fig. 5h, Step 3).With this stepwise approach, we determine the in vivo pathological molecular dynamics in urticaria patients via their eruption geometry, by extracting the geometric features of skin eruptions from patients with urticaria and linking these features to the parameter values of the CSU mathematical model.
In the first step, we successfully developed a hierarchical urticaria mathematical model in which all the model equations were verified by in vitro experiments.To our knowledge, this is the first study wherein the pathophysiological steps involved in the development of CSU were structurally proven in silico.The first conceptual mathematical model of urticaria was proposed in our recent work 7 .In the previous model, we proposed that the inhibitory regulation of histamine release from mast cells plays an important role in creating multifarious patterns of the wheal geometry, although this finding has not been experimentally verified.In this study, the inhibition dynamics of adenosine were modeled based on the experimental data 16 and proved that the inhibitory regulation in mast cells is likely to be involved in the classification of eruption types with interactions in other physiological networks (Fig. 4g, Fig. S7).Moreover, our CSU model proposes that the hierarchical structure of urticaria development as well as the key players, such as TF and histamine released from mast cells, are critical in generating diverse eruption patterns and determining the eruption types.
In the process of developing the mathematical model, we also found that the intracellular gap in endothelial cells of the blood vessels is formed at a threshold level that depends on TF concentration.This finding suggests that wheals may develop only when stimuli are presented at a level above a certain threshold.
Thus, the gap formation dynamics based on the TF switch may affect not only the formation of eruption patterns but also the onset of wheal formation.
Using our model, we also identified five eruption patterns and created EGe criteria for the classification of urticaria eruptions in real patients.Application of the EGe criteria confirmed that patients with CSU could be classified into one of the five types, and the types included in the area pattern were observed more frequently than those in the boundary pattern in the patients.Authors experience boundary patterns specifically in patients with spontaneous urticaria, whereas area patterns in both spontaneous and inducible types of urticaria (data not shown).Therefore, area patterning may reflect the ubiquitous molecular mechanisms common to patients and various types of urticaria, thereby resulting in a high patterning frequency.Taken together, these results indicate that our EGe criteria could be a tool of practical approach for clinicians to apply in determining eruption patterns of patients in the future.In the final step of our research strategy, we inferred the pathophysiological networks that are critically involved in each feature of the eruption geometry.We found that the balance of effects by TF and the histamine released from mast cells play a critical role in determining the types of eruption patterns.Our result suggests that the level of TF expression may be involved in an underlying mechanism mediating the creation of notable differences in eruption shapes such as boundary or area patterns.Furthermore, the kinetics of histamine released from mast cells may also play a role in determining the specific type of eruption in different patients.This implies that even suppressing the action of a specific network using a drug may not always exhibit the same effectiveness among patients exhibiting different patterns of eruptions.In fact, antihistamines are effective for 70-80% of CSU patients 21 , but certain patients refractory to antihistamines respond well to anticoagulants and protease inhibitors [22][23][24] .Thus, our approach, which proposes a pathophysiological in vivo network considering both topological and kinetical reactions in the skin, is expected to suggest a patient-specific treatment based on the geometry of skin eruptions (Fig. 5h, red arrows).
In the model assumed in this study, histamine plays an important role in a cascade activating mast cells.Our previous studies suggest that histamine may play both in an early step in inducing tissue factor expression on vascular endothelial cells and a late step in inducing plasma extravasation via the increase of vascular permeability (Fig. 1) 25 .However, its role in the cascade, especially in the early step, may be substituted to some extent by other factors, such as LPS and cytokines.Therefore, the inhibition of histamine by antihistamines in an early step of the cascade may or may not fully prevent the process of consequential steps involving the coagulation pathway, depending on the degrees of the involvement of histamine in this step in individual patients.In this study, we also revealed the association between high TF activation and the area pattern rather than boundary ones by theoretical analysis.Such an association may be experimentally confirmed in future studies by examining levels of coagulation biomarkers such as PF1 + 2 and D-dimer 26,27 .Thus, the relation of wheal pattern to the biomarkers of coagulation as well as autoimmunity activating mast cells and basophils needs to be evaluated in the future.
Urticaria involves both temporal and spatial dynamics, and the eruptions show highly dynamic transitions from the phases of eruption emergence (onset), development, and disappearance.The distinction between annular and broken-annular patterns, and that between circular and geographic patterns may be sometimes difficult in certain clinical practices.However, the five wheals' patterns are clearly distinguished in the course of mathematical topology (Supplementary Table S2).The wheal expands of annular patterns keep the connected boundary, whereas those of broken-annular develop with disconnected boundaries (Fig. 3).Such difference should be more apparent in kinetics information rather than a snapshot image of individual wheals.This feature is also clearly distinguishable in mathematical topology.Similarly, the features between circular and geographic patterns will be more clearly distinguished in the temporal dynamics of wheal fusion.In in silico case, we characterized the circular pattern of wheal by slow speed of expansion and sparse distribution of the emerging points.Thus, time-lapse views of wheals are expected to be obtained for analysis and endorse the model proposed in this study in future studies.
This study is limited to focusing on the emergence and development phases, and we did not consider the disappearance phase.In the current model, the eruptions would continue to expand if the simulation space is not limited.This suggests that the disappearance phase may be caused by a network independent of the onset phase.Understanding the overall dynamics of CSU from the onset to the disappearance phase will be crucial to understanding the exact pathological state of patients based on the eruption state.To achieve this, the time-related historical data of patients will need to be acquired and added to our current mathematical model.Furthermore, we constructed a mathematical model based on in vitro experimental data and succeeded in explaining the characteristics of the eruption geometry.However, a more accurate clinical and biological validation using another independent data set is needed.
Despite these limitations, the findings of this study suggest that the skin eruption geometry is linked to the dynamics of in vivo molecular networks involved in the pathogenesis of CSU.The inference of the detailed parameter values in our mathematical system for each patient may be applied to decide patient-specific treatments.Moreover, the mathematical modeling approach used in this study may be applicable not only to urticaria but also to a variety of skin diseases characterized by the formation of geometric morphology, such as psoriasis, tinea corporis, erythema multiforme, autoimmune bullous disease, erythema gyratum repens.Thus, our multidisciplinary approach linking the skin eruption geometry and in vivo, dynamics could help develop a novel strategy for dermatological research to manage intractable skin diseases in an integrated manner.
Fig. 5 Linking key features of wheal geometry to related networks.a-e Results of the first index determined using extended Fourier amplitude sensitivity test (eFAST) sensitive analysis of each feature.The first-order index of the dummy parameter was calculated to confirm that index values in panels (a-e) yielded reliable data to determine the significance of model parameters 28 .For model parameters with an index less than or equal to that of the dummy parameter, the parameters are considered not significantly different from zero.Networks with indexes lower than (negligible effect) and double the height of the dummy parameter (high effect) are indicated in blue and red, respectively.f Summary table of sensitivity analysis.Parameter sets related to basophil, TF, and mast cell dynamics are shaded by the same color.g Summary of analysis results of the pattern types and main network.The table at the bottom shows relative levels of tissue factor (TF) expression (green shaded region in the upper panel) and histamine released from mast cells (blue shaded region in the upper panel) for each eruption pattern.(h) Summary of workflow and clinical application.Black arrows show the study workflow and red arrows represent the general clinical applications of our research, leading to potential patient-specific treatments for intractable skin diseases.
SB ¼Number of doctors who selected the type in EGe criteria Total number of doctors who participated in the classification :

Fig. 1
Fig. 1 In vivo pathophysiological network of chronic spontaneous urticaria (CSU).a Schematic diagram of CSU development predicted using in vitro experiments.Number labels correspond with those in (b).Black arrows and green symbols indicate positive (activation) and negative (suppression/ inhibition) regulation, respectively.The order of pathophysiological progression is sequentially numbered.b Summary of the network based on the description in (a).Networks that do not appear explicitly in the mathematical model are depicted by gray dotted arrows.In the mathematical model, arrows (2 a ) and (2 b ) were simplified to solid arrows (2).Networks (2 a ) and (5) were integrated into the modeling of the gap-forming network (3).LPS lipopolysaccharide, H1R histamine H1 receptor, TLR4 toll-like receptor 4, PAR-1 protease activated receptor-1.

Fig. 2
Fig. 2 in vitro experiments for key networks and estimation of model functions.a, b Inhibition ratio of histamine release from human skin mast cells (hsMCs)/human peripheral basophils in response to adenosine, calculated from data of 16 .c Inhibitory functions estimated from data of a, b with the data conversion of adenosine to histamine (Supplementary Method 2).Black and white circular points are the experimental data of human peripheral basophils and human skin mast cells, respectively.The blue and red lines are the fitted graph for each case.d Inhibition ratio of tissue factor (TF) mRNA from HUVECs induced by adenosine, calculated from experimental data of Yanase, Morioke, et al., 2018.e Inhibitory function estimated from the data of (d).f Relative expression of TF mRNA from HUVECs induced by co-stimulation with histamine and LPS.Each color corresponds to the TF concentration in the right panel.g Activation by TF estimated from data of (f).κ is given in equation (s3) of Supplementary Method 3 and κ is a scaling parameter that converts LPS concentration to histamine concentration (details are presented in Supplementary Methods 2 and 3).h TF concentration-dependent gap formation.Temporal variation in cell index (CI) at different TF concentrations.Endothelial cells were stimulated with TF at indicated concentrations in Dulbecco's modified Eagle's Medium (DMEM)/F12 with 2% bovine serum albumin (BSA) and 1% factor VIII-deficient human plasma without any supplements and fetal calf serum (FCS).i Relative scale of CI, which was scaled by data obtained at 100 ng/mL each time.j Estimated switch function, χ sw TF ½ ð Þ; obtained by quantitative analysis of variation in CI at different TF concentrations.In a, b, d, f, and i, the black point indicates the mean, and the error bars are standard deviations.The blue circles in a, b, d, and f are the data points.

Fig. 3
Fig. 3 Five wheal patterns identified using a mathematical model.a-e Evolution of each wheal pattern.f Comparison of wheal pattern using a mathematical model at a time point with that of urticaria patients.Details of parameter values are shown in Table S1.g The panel shows eruption patterns by changing parameter spaces: activation of histamine release (γ M ), decay rate of histamine (μ M ), and inhibition rate of histamine (α M ) from mast cells.Yellow regions indicate the parameter region in which annular/broken-annular patterns appear.Green regions indicate the parameter region in circular/ geographic patterns.Gray regions indicate the parameter region in dot patterns.h A schematic diagram of the eruption pattern based on three parameters.

Fig. 4
Fig. 4 Classification of eruption patterns by the eruption geometry criteria (EGe criteria) for CSU patients.a Ratio of patient samples in which six dermatologists selected one of the pattern types defined by the EGe criteria.b Ratio of patient samples selected as not classifiable by the EGe criteria.Data were calculated based on reliability SB.Biologically independent patient samples (n = 105) were analyzed in (a) and (b).c Frequency of patient eruption types.Data (n = 35 biologically independent patient samples) were selected among n = 105 biologically independent patient samples with reliability SS > 83%.

Table 1
Criteria for classification of eruption geometry (EGe criteria).

Table 2
Summary table of eruption geometries based on key features extracted from EGe criteria.The table below provides a simple means of classifying eruption geometry based on five key features.