Locating Nidi for High-Frequency Chest Wall Oscillation Smart Therapy via Acoustic Imaging of Lung Airways as a Spatial Network

High-frequency chest wall oscillation (HFCWO) therapy is one of the techniques to facilitate the draining of a patient’s lung secretion in pathological situations, and smart therapy with HFCWO devices equipped with multiple actuators can be achieved via locating nidi in the lung. In this paper, through developing a novel acoustic lung spatial model and utilizing acoustic imaging simulation, a new and effective method for assessing lung function with acoustic imaging is presented, which links acoustic lung images with pathologic changes. The structural similarity between the acoustic reference image based on actual lung sound and our model acoustic image based on the airway impedance was achieved by an index of 0.8987, with 1 as the exact score. Simulation studies based on the model are used to analyze the practicality and the extreme design of the acoustic imaging system on the resolution of the located nidus. For instance, a practical system design with sensor numbers between 4 and 35 may recognize a lower resolution nidus length of 73 mm to a better resolution nidus length of 22 mm. On the other hand, an extreme system design with more than 1000 sensors can recognize greater nidus resolution at under 10 mm. Additionally, this research may be utilized to offer recommendations for acoustic imaging system design and assess the number of sensors and sensing diameter in current acoustic imaging systems. Furthermore, the geographic detection of nidus length allows for analyzing of HFCWO therapy results.


I. INTRODUCTION
Chronic inflammation, cystic fibrosis, and some respiratory viral diseases cause mucous discharge to thicken.High-Frequency Chest Wall Oscillation (HFCWO) therapy is a common airway clearance technique for patients with thick mucus and low mucociliary clearance (MCC) efficiency.HFCWO devices are defined as small oscillations of mechanical parts at relatively high frequencies (5-20 Hz) applied onto the patient's thorax for respiratory therapy.Traditional The associate editor coordinating the review of this manuscript and approving it for publication was Mohammad Zia Ur Rahman .
HFCWO devices, such as the Vest 105 by Hillrom [1], use an air-filled garment enclosing the patient's chest to generate motion similar to MCC.The parameter setting and operation are purely empirical according to user experience.Modern HFCWO devices such as the Monarch [2], the AffloVest [3], and the RespIn 11 [4] were equipped with multiple electromagnetic/pneumatic actuators that can be controlled individually, enabling a smart therapy that targets the nidus locations for an optimal therapeutic process.Therefore, knowledge of nidus location in the airway is critical.
This study presents studies on acoustic imaging to locate nidi to allow inference on the efficiency of HFCWO physiotherapy by respiratory remodeling and acoustic imaging sensor array design simulations.To the best of our knowledge, locating nidi through the two-dimensional (2D) acoustic lung model and the resulting acoustic lung imaging have yet to be performed.Moreover, the acoustic imaging system setups are typically empirical, potentially leading to unoptimized nidus detection.Hence, the key contributions are: 1) Proposing a realistic 2D acoustic lung model incorporating spatial location to simulate airway obstruction and to design and optimize acoustic sensor array measurements quantitatively [5], [6], [7], [8].2) Applying the resulting acoustic image from the proposed 2D airway model to theoretical acoustic sensor array design by considering the sensor distribution, sensor sensitivity area, and sensor number.
First, by predetermining the acoustic sensor sensing area, this research illustrates the relationship between the severity of the airway obstruction and mean acoustic image intensity through the thickening of the airway wall thickness (AWT).A good agreement was found between a reference obstructed airway created from lung sound data and acoustic imaging from our model, with a structural similarity (SSIM) index of 0.8987, with 1 denoting an identical image.Next, different sensor sensing areas are employed to correlate the observed nidus length with the sensor numbers.About 26,000 sensors are required to identify a resolution of 4.35 mm minimal nidus length with a 10 mm sensor sensing diameter.Comparatively, a 50 mm sensor sensing diameter may identify a roughly 73 mm minimal nidus length resolution with only about 4 sensors.The findings support the theory that better image resolution derives from increased sensor numbers.In addition, the required sensor numbers and sensing sensitivity can be used as a baseline consideration in the acoustic imaging system design.Additionally, a guideline for designing HFCWO devices and assessing the HFCWO therapy efficacy on the patient for a smarter process through therapy feedback from identified nidus length can potentially be provided by understanding how sensor array and sensing sensitivity affects lung health assessment with the resolution of detected nidus and optimizing the sensor array.
This paper is organized as follows.An incisive review of the airway modeling and the acoustic sensor array design are presented in Section II.The modeling of airways and generation of the acoustic imaging are described in Section III.Model verification by comparing healthy lungs and the lungs with asthma and chronic obstructive pulmonary disease (COPD) symptoms are demonstrated in Section IV.The simulation studies on locating nidi, sensor distribution, and image resolution are presented in Section V, followed by general discussions in Section VI.Lastly, the conclusion and future work are given in Section VII.

II. LITERATURE REVIEW
To realize the HFCWO smart therapy, locating nidi is critical, while one of the direct ways to access nidus location is to present on an image.Chest X-rays, CT, and magnetic resonance imaging (MRI) are the usual imaging techniques to visualize the airways and lung pathology.However, these approaches are not ideal due to their ionizing radiation effects and the 'patientto-equipment' approach [9].Unlike chest X-rays, CT, and MRI, electrical impedance tomography (EIT) [10] is an 'equipment-to-patient' approach and uses nonionizing radiation technology that provides alternatives to monitor airways.However, EIT usually provides transverse plane images instead of the required frontal plane images (see Fig. 1(b)) for the actuator selection or adjustment (see Fig. 1(a)), making it challenging to apply to HFCWO therapy.
In the quantitative forms of lung sound presentation, Kompis et al. [11] developed an acoustic imaging technique that uses simultaneous multimicrophone recordings to assess spatial information.Another technique for converting the acoustic signal to an image is Vibration Response Imaging (VRI) [5].VRI reflects the dynamic changes in the lung by imaging that utilizes the vibration energy created during breathing.By presenting localized information on breath sounds between different lung sites, the visual representation improves the clinical value [5].Acoustic imaging and lung disorders, such as smoking index and the accumulation of extra fluid between layers of the pleura outside the lungs, have a positive quantitative data correlation [5].Computing from the impedance or the resistivity in the lung or the airway through respiratory remodeling as an indicator for lung function assessment is required as an initial step.
Airway obstruction or the thickening of airway wall occur in chronic respiratory illness, alter the production and transmission of lung sound spectrally and regionally.Asthma and COPD patients with frequent mucus production in their airways tend to have thicker airway walls than those without, regardless of the severity of breathlessness, and have shown significantly different morphologic airway findings compared to healthy individuals [12].The change can be measured quantitatively in the lung sound transmission and provide critical information on the disease severity and location of the airway obstruction [12], [13], [14], [15], [16], [17].Spatially distributed airway tree models have been developed to decipher the relationship between bronchi lengths, branching angles, and airway diameters [18].In the development, Murray's law [19] defined that the relationship between airway bifurcation is fixed, with branch lengths based on a length-to-diameter ratio.Weibel symmetric and Horsfield asymmetric models are the most used conducting airway models [18].With the advancement of medical imaging techniques, deterministic parameterized bronchial tree generation algorithms were extracted directly from computed tomography (CT), thus constituting the core of patient-specific modeling [18].The recent works in this area are summarized in [20].However, those models developed so far are typically simplified to a one-dimensional system of equations to investigate the relationships between healthy and unhealthy respiratory system cycles, such as frequency response, flow rate, resistance, volume, and diagnosis accuracy [13], [18], [21].
Although positive correlation can be identified through acoustic imaging and lung disorders, chronic respiratory diseases, such as asthma and COPD, have not been correlated positively [5], [22].Moreover, the correlation between acoustic sensor placement and sensitivity was not investigated, and the position of the sensors was typically empirical [5], [11], [23].The summary of the key points and the research gap identified from the concise literature review is presented in Table 1.

III. MODELING OF AIRWAYS AND ACOUSTIC IMAGING
This paper developed a model for acoustic imaging with the following features to improve the investigation of locating airway obstruction, as each patient has a unique set of airway dimensions and structures: 1) The ability to modify the airway input parameters that influence the model's output, such as the wall thickness, length, and diameter, where the patientcentric assessment technique is made possible.
2) The airway model outputs intuitive spatial-based 2D imaging to show airway obstruction in the lung caused by respiratory conditions such as COPD and asthma (Section IV).
3) The resolution of the lung image was intended mainly for the assessment and location of the obstruction in the airways due to the limited sensor numbers and HFCWO actuators that can fit onto the patient's posterior chest area (Section V).
Drawing inspiration from [11], [13], [15], [16], [20], and [21], the respiratory system is represented as a bifurcating tree network with the linked node of the bifurcating segment and integrated spatial position (x, y) on the airway plane, where the airway plane refers to the three dimensional (3D) airway network space that is projected onto.After that, the network is converted into an electrical network with lumped characteristics and presented as an assessment of the acoustic lung image.In the model development, the following notations are used.R denotes the set of all real numbers.R m×n is the set of all real (m × n) matrices.C denotes the set of all complex numbers.C m×n is the set of all complex (m × n) matrices.Z(ω) is the set of all sinusoidal variables with angular frequency ω.
The construction of respiratory airway modeling on a single node of the bifurcating airway impedance and the respiratory airway modeling parameter is presented in Sections III-A and III-B, respectively.Next, the conversion of the airway impedance into acoustic imaging is presented in Sections III-C.

A. MODELING RESPIRATORY AIRWAY
Each 3D network segment is initially projected toward a 2D plane and given a coordinate for its position (x, y).The respiratory system is thus depicted as a bifurcating tree network, with the joined node of the bifurcating segment at layer k and position (x, y) being indexed by (x, y, k) on the plane illustrated in Fig. 2(a).Through a recursion index of (k), the k-th layer segment splits into asymmetrical airways of layers (k + 1) and (k + 1 + (k)) [15].The airway is then represented as a network of bifurcating cylinders which can be modelled as a transmission line with distributed parameters and further translated into an electrical network with lumped parameters, as shown in Fig. 2(b).The airway network is then resolved by the acoustic pressure at each segment induced by the pressure distribution from bronchi breathing and the airway network [24], [25].Merging the acoustic power over a predetermined period of time during each breathing cycle, a plane image is generated by the projected network as a subset of the acoustic lung image Q(x, y) ∈ R m×n (discussed in Section III-C).
Since the longitudinal motion of the airway is typically negligible in comparison to the acoustic signal, the acoustical impedance Z (ω) and acoustical admittance Y (ω) averaged where l is the axial coordinate.When the patient breathes periodically, the airway can be regarded as a steady-state system with each segment as a short nonrigid transmission line tube with unit-length parameters equivalent acoustic resistance R 0 , inductance L 0 , capacitance C 0 , and conductance G 0 [12], described by (2), where P 1 ∈ Z and F 1 ∈ Z are the input pressure and input flowrate and P 2 ∈ Z and F 2 ∈ Z are the output pressure and output flowrate, respectively.The propagation coefficient γ The transmission line tube with distributed parameters can be equivalent to a network in Fig. 2(b) with lumped parameters of segment impedance Z g ∈ C and segment admittance Hence, the entire network of airways can be represented as an electrical network made up of a layered bifurcating tree of impedance connected to the ground through an admittance at each bifurcating node, as illustrated in Fig. 3(a).The air pressure and airflow rate are comparable to electrical potential and current, respectively, when the respiratory airways are analyzed as an electrical network [13], [21], [26], [27], [28], [29].The k-th layer's impedance and admittance can be presented in (5), The network of airways is constructed with n nodes indexed with encircled numbers, b branches denoted with underlined numbers, the k-th layer as subscript, and a sinusoidal voltage source with amplitude P s and angular frequency ω in series of a small impedance Z s0 applied at the input layer 0 to represent the fundamental component of the periodical patient breath, as presented in Fig. 3(a).An incidence matrix A will be used to evaluate and simulate an acoustic network encompassing resistive and capacitive elements scattered over multiple interacting layers and acquiring a descriptor representation of the network, as demonstrated in Fig. 2, Fig. 3, and (1)-( 5).Thus, we have the following annotation shown in (6) from the theory of network topology [24], [25].
where A, Y, Y b , V s , V, V n , I s , and I, are reduced incidence matrix, branch admittance matrix, branch admittance vector, branch voltage source vector, branch voltage vector, node voltage vector, branch current source vector, and node current vector, respectively.A standard branch in a linear network is shown in Fig. 3(b), and the node analysis is given in (7), where Kirchhoff's voltage law and Kirchhoff's current law serve as the first and second requirements in (7), respectively, with the third requirement deriving from the standard branch law, and (8) can be obtained from the node analysis in (7).
The node voltage V n is the remaining unknown variable from (8).Assuming node admittance Y n ∈ C (n−1)×(n−1) is a nonsingular and symmetric square matrix, and J s ∈ Z (n−1) is the node source-current vector as shown in ( 9), the node voltage V n can be resolved in (10), From the graph in Fig. 3(a), assuming b = 3 × 2 n , and I s = 0 in ( 6)-( 9), the reduced incidence matrix A and branch admittance matrix Y can be denoted as follows: Table 2 shows the incident matrix A, branch admittance vector Y b , and branch voltage source vector V s of the first four network layers in Fig. 3.The reduced incidence matrix A is the resulting network matrix without the row of node G in Table 2.
Given that the patient's breath pressure is sinusoidal, every joint pressure can be resolved by the network analysis method as long as parameters Z k and Y k are known.From (5), this needs to find the segment parameters Z g and Y g .

B. PARAMETERS OF RESPIRATORY AIRWAY MODEL
The airway wall was modeled using the complex Young's modulus and material density to replicate the acoustic structural interaction accurately [13], [15], [21], [30], where the material parameters of the respiratory system are given in Table 3.The airway segments' thickness, cartilage, and soft tissue fractions were determined by referring to the data reported in [15] and identifying the closest Horsfield order segment.Thus, the segment in the k-th layer has the material parameters in (12), where A s (k), α(k), ρ g , η g , C g , K g are denoted as the cross-sectional area of an airway segment, internal airway radius, air density, viscosity, air specific heat, and thermal conductivity, respectively [14], [27], [31].F v (k, ω) and F t (k, ω) account for the sound attenuation by air viscosity and sound attenuation by thermal dissipation, computed with series expansion with J 0 (z v ), J 1 (z v ) and J 0 (z t ), J 1 (z t ) being Bessel functions of 0-th and 1-st orders [31].Z w (k, ω) represents the wall impedance, which is computed from a series of resistance R i (k), inductance L i (k), and conductance G i (k) of the acoustic transmission line and Young's modulus E i , where the subscript i is replaced by either c for the cartilage or by s for the soft tissue, respectively.

C. ACOUSTIC IMAGE GENERATION
Most of the previous works investigate the variable physical frequency characteristics [13], [21], [26], [27], [28], [29], and no spatial information is associated with the nodes.In this study, the spatial location (x, y) was integrated into each node to transform the airway network into a spatial network and generate the resulting acoustic image.The acoustic image can be initiated once the node voltage V n , which is analogous to the acoustic pressure P n distribution within the airways [13], [21], [26], [27], [28], [29] is obtained.The sound pressure in dB within the airways is computed as, P = 20 log 10 (P n /P 0 ), (13) where P 0 = 20 µPa is the reference sound pressure.
The sound pressure generated from the lumped electrical network resulting from the transformation of the respiratory modeling, as presented in Fig. 2 and Fig. 3, can be captured with an array of acoustic sensors (see Fig. 4), such as digital stethoscope or micro-electromechanical systems (MEMS) microphone [5], [11], [23].An interpolation function can be utilized to compute the sound pressure between each sensor [23].
The airway pressure at each sensor location is computed by accumulating the captured signals over a given time interval t from t 1 to t k and averaging the signals at all bifurcating airway nodes within the sensing area enclosed by the horizontal and vertical lines as the individual area boundary as shown in Fig. 4, where N s is the total number of airway nodes within the sensing area.The network of the acoustic lung image Q(x, y, t 1 , t k ) is then, The sound intensity outside of the sensor position in Fig. 4 is estimated by interpolation.From the observation in ( 13)-( 15), the acoustic lung image Q( P, h) is defined as the 2D acoustic image which comprise acoustic signal P(x, y, t 1 , t k ) in ( 13) and ( 14) with interpolation polynomial h(x, y).A high spatial resolution is required; hence, Hermite interpolation was applied to the acoustic signal P for projecting acoustic lung imaging [23].From the study in [23], Hermite interpolation has been proven to be a better performance in presenting accurate lung sound intensity as compared to other established interpolation functions, such as linear, cube spline, Lagrange and nearest neighbor method.Refer to [23] for the Hermite interpolation function in-depth analysis, computation and application on acoustic lung imaging.
Each acoustic image pixel is normalized, and the output obtained from the pressure sound signal is then presented as an acoustic image with the highest, lowest, and in-between values are determined as maroon, white, and grey.

IV. MODEL VERIFICATION BY PATHOLOGY EXAMPLES
Model verification and the potential to assess the severity of airway obstruction through regional pathology with a predetermined sensor number and sensor sensing area are demonstrated in this section.Additionally, due to the vast range of lung sound frequencies documented in the literature, 400 Hz was chosen as the frequency to convey the results in this paper for the relevancy to respiratory sounds and to keep it straightforward [21], [32].
A reference image was produced from a COPD patient's lung sound signal that was selected from a respiratory database [33].A four-by-six array of sensors, as illustrated in Fig. 4, where the sensors are considered to be equally dispersed within a 50 mm distance [5], [11], [23], and the acoustic response is the average intensity value within the sensing region.A 2D plane acoustic lung image can be produced with ( 13)-( 15) and the known sensor and spatial position information (x-and y-axis) as shown in Fig. 2 and Fig. 4. The light-colored (white) area is the color for the minimal or no pressure data area, which represent the airway's high airflow resistance, whereas dark-colored (maroon) area is used to indicate high data locations where the airflow resistance in the airway is the least.Additionally, the in-between data area, where the airway has airflow resistance, is represented by light gray colors.
In the following, the assumption for the model simulation of pathology through AWT remodeling and the quantitative model performance are presented in Section IV-A and Section IV-B, respectively.The results and discussion are described in Section IV-C.

A. PATHOLOGY SIMULATION
Airway remodeling was performed by altering the AWT to simulate airway obstruction [34], [35], [36], [37], [38], [39].As shown in Fig. 2, the total wall thickness of each airway segment H w = D o − D i , where D i and D o are the inner and outer diameters, respectively.The inner airway diameter D i and total wall thickness H w were measured and compared from patients with illnesses, such as asthma and COPD, using computed tomography in relation to the severity (mild, moderate, severe) of the illness [34], [35], [36], [38], [40].The studies in [34], [35], [36], [38], and [40] have revealed a range for the mean airway wall area percentage (WA%) increment of 3%-40%, with 0%-3% for controls, 4%-10% for mild conditions, 11%-30% for moderate conditions and more than 30% for severe conditions.The studies on airway wall thickness and the increment of the airway wall area are summarized and presented in Table 4.The airway wall area (WA) and WA% can be calculated as [35], where A 0 = π(D o /2) 2 and A l = π(D i /2) 2 can be computed as the airway area and the luminal area, respectively.

B. PERFORMANCE ASSESSMENT
The mean acoustic image intensity (dB) in ( 13)-( 15) can be utilized as an indicator for the assessment outcome on the severity of airway obstructions [5], [6], [7].The increment (factor) of AWT was implemented to standardize the findings in this study, as mixed airway obstruction results can be identified from the literature, such as the increment of WA% or values of AWT [34], [35], [36], [38], [40].
For instance, the AWT must increase by a mean factor of 2.34, as shown in (16), for the mean WA% to increase by approximately 11%, from 67% healthy lung to 78% respiratory illness lung [34], [35], [36], [38], [40].Finally, in terms of the severity of respiratory diseases, the internal airway area between asthma and COPD was essentially the same [34].Therefore, no differentiation between COPD and asthma is made in this study.
The pixels in natural image signals are heavily dependent on one another, especially when the pixels are close together.These dependencies include important details about how the elements in the visual scene are arranged.The SSIM index [41] is a straightforward approach for comparing the reference and distorted signal structures.Additionally, SSIM indexing provides quality assessment from the perspective of image generation, particularly for components of medical images in pixel intensities [42].The SSIM quality assessment index is based on the computation of three terms, namely the brightness term, the contrast term, and the structure term, as illustrated in (17), where µ Qr and µ Q are the local means, σ Qr and σ Q , are the standard deviations, σ QrQ cross-covariance, and C 1 and C 2 are the constants for reference image Q r and captured image Q.For detailed derivation and computation, see [41].

C. MODEL VALIDATION
The acoustic lung imaging Q projected from lung signals is computed from the lung signal intensity P at each sensor location in a coordinate plane over a known time t interval, as shown in ( 14) and ( 15).The lung signal intensity is determined as highest (maroon), lowest (white) or in between values (grey).The acoustic signal is normalized, and the output obtained from the intensity of the sound signal is then displayed as an acoustic image.The overview for the model validation in this study is presented in Fig. 5, where the increment factor for the airway wall thickness has been discussed earlier in Section IV-B, the computation of the acoustic signal can be identified from ( 1)-( 12), and the translation of the computed acoustic signal to acoustic image can be inferred from ( 13)-( 15).An unaltered airway was utilized as a control in this model validation, where the airway material properties and parameters have been introduced earlier in Section III-B.The spatial resolution of the lung geometry in this model validation is 44 pixels for every 10 millimeters.Fig. 6 displays acoustic images of a healthy lung (control) and varying respiratory illness severity obtained by adjusting the AWT in Section III and ( 13)-( 15).An outline is used in Fig. 6 to identify better the effect of AWT on the overall (global) lung image intensity.Additionally, Fig. 7 displays the relationship between the average image intensity and the global AWT increment.
Fig. 6 and Fig. 7 demonstrated the relationship between acoustic lung images of healthy and ill conditions.In contrast to ill conditions, such as mild, moderate, and severe conditions, a healthy lung presents the darkest lung image (high acoustic intensity value) due to the lowest impedancesmallest resistance in the airway, from Fig. 6 and Fig. 7 and as observed from ( 2)- (13).Moreover, the airflow and the mean image intensity both reduced with the thickening of AWT can be observed in Fig. 6 and Fig. 7.Although observable qualitative changes can be seen with the AWT  increasing by a factor of more than 1.70 in Fig. 6, the mean image intensity in Fig. 7 can reveal the state of the lungs' condition.Furthermore, the positive correlation between the lung impedance from ( 2)-( 13) and the results in Fig. 6 and Fig. 7 presented a certain level of similarity compared to the literature [5], [6], [7], [15], [34], [35], [40], e.g., the global intensity distribution impacting the lung and the airway closer to the trachea (Fig. 2) is often larger and tends to be the last impacted region by the thickening in AWT.
After the global thickening in AWT and the consequences (severity) on lung function have been demonstrated, the next validation task is the regional increase in AWT.Fig. 8 contrasts our model acoustic image with the obstructed reference lung image, which was created using the lung sound signals extracted from a respiratory database [33], and converted into an acoustic image.The obstructed airways are situated along the posterior right middle scapular line (area B2), and the posterior right lower scapular line (area C3), as shown in Fig. 4. The region of the obstructed airway can be located in our model's acoustic image presented in Fig. 8.The similarity between the acoustic reference image and the model acoustic image is highly related given that a mean SSIM index of 0.8987 was obtained, with 1 being the same as [41].

V. OPTIMAL ACOUSTIC SENSOR ARRAY DESIGN FOR AIRWAY OBSTRUCTION DETECTION
Global and regional pathology with prearranged number of sensors, e.g., an array of 4-by-6 with 50 mm uniform spacing acoustic sensors, has been validated in Section IV.The remaining task in this study is how the design of the acoustic sensor array affects the minimal detectable nidus length, e.g., the expected minimal detectable nidus if the acoustic sensor array is known or the design of acoustic sensor array for an envisioned minimal detectable nidus length.To the best of the authors' knowledge, no discussion was attempted relating to the distribution array of acoustic sensors for image assessment and the acoustic imaging resolution, as the array sensor design was typically empirical in the literature [5], [21], [23].In line with the uniform distribution design of HFCWO electromagnetic/pneumatic actuators [2], [4], and the traditional acoustic imaging system in the literature [5], [21], [23], a uniform multimicrophone distribution, vertically and horizontally, is employed in this paper.In addition, the overlapping and nonoverlapping sensor sensing sensitivity can be computed due to the influence of the sensor uniformly distributed.Hence, the effect of sensor sensing sensitivity area and the sensor number on the detection of airway obstruction is presented in Section V-A, followed by analysis and discussion in Section V-B.
By employing local first-order image statistics [43] around each pixel, the resulting obstructed airway acoustic image µ are converted into a binary image, as shown in Fig. 9.As shown in Fig. 9(c), areas with high-intensity data (healthy) are denoted by 1s, and areas with low-intensity data (obstruction) by 0s.Thus, by comparing the acoustic image pixel area η in Fig. 9(a) and the pixel area µ in Fig. 9(c), the obstruction in the airway acoustic images can be located, and the area of the missing pixel (η − µ) can then be used to calculate the obstructed area (nidus) length,

A. SENSOR SENSING SENSITIVITY AND SENSOR NUMBER
To study the effect of sensor sensitivity on the smallest observable nidus length L n , the number of sensors is initially fixed at 12, 16, 20, 25, 32, 40, 45, and 50 per lung side, comparable to the empirical acoustic image system [5], [21], [23].The selection of the sensor sensing diameters, which ranged from 10 mm to 50 mm in 10 mm increments, was made in accordance with commercially available products and published research [5], [23].Fig. 10 shows the relationship between sensor sensitivity with a predetermined number of sensors and the measured minimum nidus length.
After the effect of the different sensor detecting areas on the minimum detectable nidus length when used with a predetermined number of sensors, the next step is to evaluate how the number of sensors affects the minimal detectable nidus.Fig. 11 illustrates how the number of sensors affects the minimum observable nidus length for different sensor sensing sensitivities.

B. ANALYSIS OF THE SENSOR ARRAY DESIGN
A minimal detectable nidus length of about 68 mm is expected when using 12 sensors with a 10 mm sensor sensing diameter, as illustrated in Fig. 10.In contrast, a minimal detectable nidus length of about 20 mm is expected with 50 sensors, with sensor sensing diameter between 20 mm and 50 mm.Fig. 11 demonstrates the number of sensors and the sensor sensing diameter required in the acoustic sensor array for the envisioned minimal detectable nidus length.
According to Fig. 10, a better resolution of the detectable minimal nidus length was obtained with the increase in sensor number, and the sensor sensitivity area overlaps more when compared to fewer sensor numbers and lesser overlapping of sensor sensitivity area.When compared across all sensor  sensing diameters, the predefined sensor number showed various observed nidus lengths, as shown in Fig. 10.The results are in line with the number of sensors and the position, where higher image resolution can be identified with sensor sensing diameter (30 mm-50 mm) overlapping reducing the over-reliant on interpolation function, as compared to sensor sensing area that has lesser nonoverlapping sensor sensing diameter (10 mm-20 mm) [23].
From Fig. 11(a), a resolution observed in the detected nidus length is about 73 mm, requiring about 4 sensors, with a 50 mm sensor sensing diameter and a 0% sensor sensing overlapping area.In comparison, a high resolution identified in nidus length is about 4.35 mm, requiring about 26,000 sensors, with a 10 mm sensor sensing diameter and a 95% sensor sensing overlapping area.The observations in Fig. 10 and Fig. 11, where the resolution of nidus length detected increases with the increase in sensor numbers and corresponded with the understanding that image resolution increases with the number of sensors.The practicality in the designing of an acoustic imaging system for the location of nidus length, in terms of the number of sensors required, is demonstrated in Fig. 11(b).

VI. GENERAL DISCUSSION
The severity of respiratory diseases has been demonstrated with the mean image intensity and the thickening of AWT.The assessment of lung function through acoustic imaging, such as presenting global and regional obstructed airways, was demonstrated in Fig. 6-Fig.8.The majority of earlier studies [13], [21], [26], [27], [28], [29] focus on the changeable physical frequency features and no geographical information was correlated with nodes in the airway.Thus, the airway network in this study was converted into a spatial network by integrating the spatial position (x, y) to each node, which can produce the acoustic image for lung function assessment as shown in Fig. 4, Fig. 6, Fig. 8, and Fig. 9.All 35-airway segment layers, starting with the trachea at k = 1 and terminating at the terminal bronchiole with k = 35, were included in the calculation of the acoustical impedance.A similarity rating of about 89% was achieved between our model image and a reference image converted from lung sound signals.Minimal differences in Fig. 8 and the SSIM rating are expected as the acoustic images in Fig. 8 were generated from two different sources: our model computed acoustic impedance and the actual acoustic signal from a respiratory database [33].Only large airways, e.g., airway segment length > 2 mm, were utilized in our acoustic imaging, as small airways length ≤ 2 mm flow is laminar and silent, hence, do not produce an acoustic signal [32].Bifurcate node angles of the airway system assumed to be between 45 and 60 degrees and were drawn ideally in Fig. 2(a) so that the airway system does not overlap [21].The sound pressure computation is based on the mean sound pressure within the sensing region, as shown in Fig. 4 and (14).Hence, the SSIM rating can be improved with the additional weighted ratio between pressure in the individual airway segment and sensor sensing radius to ( 13)- (15), and an increase in the total number of airway segments in the model.
In addition, this study's objective demonstrated the respiratory model systems' capability to pinpoint the source of airway obstruction through acoustic signals, in terms of the minimal nidus length identified through the location of obstructed airways to both the acoustic sensor sensitives and the number of acoustic sensors to improve HFCWO therapy in Fig. 6-Fig.11.Although the findings in Fig. 6-Fig.11 are based on a uniform distribution of sensor location, this paper can be used as a starting point to study nonuniform sensor distribution, which may potentially result in a reduction in the number of sensors needed to achieve the same performance.Additionally, this work uses respiratory remodeling and sensor array simulation to evaluate the sensor's placement, sensitivity ranges, and the numbers for minimal nidus length detection, enabling deductions about the efficacy of HFCWO physiotherapy with the detected nidi.This paper can also be used to assess an existing acoustic array system and provide direction for the development of acoustic imaging systems, particularly in imaging systems that employ a multi-acoustic sensor array.Therefore, by comprehending how sensor array and sensing sensitivity affect lung health assessment with the resolution of detected nidus and optimizing the sensor array, a guideline for designing HFCWO devices and assessing the HFCWO therapy efficacy on the patient for a smarter process through therapy feedback from identified nidus length may be provided.A comparison to summarize the key advantages and limitations of the previous work and this study is presented in Table 5.

A. DESIGN CONSIDERATION OF IMAGING HARDWARE SYSTEM
Two of the many deciding considerations in creating the acoustic imaging systems in this study can be sensor type and sensor costs.Different transduction techniques, such as condenser (MEMS microphones) and piezoelectric (digital stethoscope) transduction, can be used to record the acoustic images derived from acoustic lung signals.Piezoelectric sensors were often not mechanically durable and required hard, specialized contacts with the patient's skin, such as gels and vacuum seals [5], [44].Due to their repeatable frequency response and high SNR, MEMS microphones are frequently employed to acquire lung sound signals and indirectly provide excellent acoustic imaging [45], [46], [47], [48].Additionally, flexible multisensor arrays, such as MEMS microphone arrays, are perfect for delivering a 2D visualization assessment of the lungs in contrast to a single sensor, such as a digital stethoscope, which can only provide one region of data at a time [45], [46], [47], [48].
MEMS microphones are also small, light, and inexpensive, costing only a few dollars, around USD 4, as opposed to a digital stethoscope, which may run between USD 300 and USD 500 [45], [46], [47].In addition, MEMS microphones can be redesigned to accommodate various sensor sensing diameter requirements, such as 10 mm, 20 mm, or 50 mm, while the sensor sensing diameter is designed to partially integrate over the fixed surface area (50 mm) of the stethoscope head [48].
For the same detected minimal nidus length, several sensor numbers and sensitivity combinations can be perceived in Fig. 11(a).A minimal nidus length of around 50 mm that can be detected, for instance, can be achieved using 6 pieces of 50 mm sensor sensing diameter or 16 pieces of 10 mm sensor sensing diameter.Given that one MEMS microphone can cover a 10 mm sensing diameter and five MEMS microphones can cover a 50 mm sensing diameter [48], using a 10 mm sensor sensing diameter may cost the customer roughly USD 64 as opposed to USD 120 with a 50 mm sensor sensing diameter.Similarly, a minimal nidus length of around 30 mm that can be detected, for instance, can be achieved using 20 pieces of 50 mm sensor sensing diameter or 48 pieces of 10 mm sensor sensing diameter.Using a 10 mm sensor sensing diameter may cost the customer roughly USD 192 as opposed to USD 400 with a 50 mm sensor sensing diameter [48].In terms of the standard MEMS microphone physical size and the adult chest area, a maximum of roughly 1000 pieces of MEMS microphone with a 10 mm sensor sensing diameter can be fitted without physical devices overlapping onto the chest region [45], [46], [47], [48], [49].We anticipate that as sensor technology advances in terms of the physical size, allowing the number of sensors to multiply, the resolution of the detectable nidus length can also be enhanced, as depicted in Fig. 11(a).Since the lung assessment imaging gold standard, such as chest X-ray, has a high operational cost (> USD 5000) and radiation exposure (health hazard), which indirectly leads to the unsuitability in frequent assessment, the detection of obstructed airways by acoustic imaging represents a crucial clinical need [44].In the literature [5], [17], [44], an array of microphones was employed to produce an acoustic image that was comparable to a chest X-ray in terms of sensitivity, specificity, and intra-and inter-rater agreement.Besides, the resolution of the acoustic lung image in our paper was primarily designed to enable frequent nidus detection by simple 2D image viewing and frequent evaluation of the efficacy of HFCWO therapy.
The computerized respiratory sound analysis (CORSA) recommendations for sensor properties to detect human pulmonary sounds can be used to guide the choice of the MEMS specification [44], [50].Other MEMS have been used to record breathing patterns and respiratory rate, a feature that can also offer a thorough analysis of lung signals.Examples of these MEMS include MEMS accelerometers [51], [52], MEMS piezoelectric resonant microphones [53], and MEMS strain gauges [54].As this study focuses on proposing a realistic 2D acoustic lung model incorporating spatial location to simulate airway obstruction and to design and optimize acoustic array measurements quantitatively by applying generic acoustic sensor array design by considering only the sensor distribution, sensor sensitivity area, and the sensor number, readers who are interested in the fabrication of the various state-of-the-art MEMS can refer to [51], [53], and [54] and the references therein for in-depth details.

B. LIMITATION
With the current study, four critical points should be considered.First, this study focused on lung acoustic signals generated from the proposed model, while the separation of heart sound signals and lung sound signals was not considered.The signals obtained were assumed to be at the patient's posterior, similar to how a doctor and clinicians perform auscultation, significantly minimizing the interference from heart signals.Likewise, the reference acoustic image translated from the actual lung signals from the respiratory database were recorded on the patient's posterior to ensure that the heart sounds would be minimal and would not significantly interfere with the lung sounds.Additionally, the frequency range for heart signals is typically below 150 Hz, while the frequencies of interest for lung signals range from 250 Hz to 1000 Hz [48], [55], [56], and 400 Hz was utilized in this work.Thus, a straightforward approach is to implement a high-pass filter to eliminate the lower heart signal frequency.Second, there will be variations in respiratory system model performance due to a range of factors such as the system network architecture: node position in the x-and y-axis location, and the physical airway model, e.g., Horsfield or Weibel airway model.The results presented in this paper are based on the respiratory model's independent abilities to optimize both the number and position of acoustic sensors for obtaining useful acoustic information, and other unsupportable combinations of acoustic sensor's position are not taken into account, such as imbalanced position, e.g., an offset position from adjacent sensors.Although breathing patterns and respiratory rate with respect to lung signals can be utilized for a more comprehensive lung function assessment other than acoustic lung sound signals, the frequency ranges for various breathing patterns and bodily movements overlap, and significant techniques to signal processing are required to isolate the signal components while restoring the important data for assessment purposes [51], [52].Third, the diameter of the obstructed lung region estimated from a circle's surface area is used to establish the length of the obstructed airway reported in Fig. 10 and Fig. 11.The airway geometry was assumed to be translated from a 3D space to a 2D plane without any intersections.To prevent outliners from determining the nidus length, a carefully selected simulated obstructed area was used.The lung size [49] of the respiratory system model shown in Fig. 2 is maintained at roughly 240 mm (height) by 100 mm (width), which is within 90% of the actual lung size.Finally, it is possible to locate the obstructed area in the simulated lung model precisely due to 1) only sensor distribution and sensor sensitivity area were considered in the simulated acoustic imaging sensor array design, and the actual sensor characteristics were excluded; 2) The model is believed to be interference-free from body movement, body temperature, ambient, and the ideal sound pressure can be captured directly through typical acoustic sensors utilized for capturing lung sound signals [45], [46], [47], [48].

VII. CONCLUSION AND FUTURE WORK
A spatial network of the respiratory system modeling is proposed in this paper, and sensor array design studies through acoustic lung imaging based on the model are conducted.The study results in a framework for the optimization of the HFCWO therapeutic technique that has shown: 1) The acoustic relationships and imaging characteristics between the sensing system and the location of nidus; and 2) How the sensor numbers and sensor sensing sensitivity affect the image dynamics at various locations within the chest area.The potential of assessing lung function with acoustic imaging has been validated through respiratory remodeling and obtained a similarity of 89% as compared to the acoustic image initiated from actual lung sound signals.This study offered design guidelines for acoustic imaging systems, or served as a performance assessment of already-in-use multimicrophone array-based acoustic imaging systems.Although there have been experimental studies on the location of the nidus, these researches in [13] and [21] concentrated on acoustic sound detection rather than acoustic imaging and did not take into account the impact of sensor sensitivity or sensor number [8], [13], [21].In order to support the conclusions in Sections III, IV, and V about respiratory system modeling and sensor array design, an experimental investigation on locating nidus using an acoustic imaging system can be carried out.Lastly, this work can be further used to compare the modeling and simulation results with the actual respiratory lung sounds that contain noise interferences.Thus, as a long-term goal of this research, it is possible to investigate the impact of nonuniformly distributed sensor configuration on nidus detection and the addition of a denoising algorithm [56], [57] to the acoustic imaging system for a practical system to precisely identify the location of the pathology produced by the airways for targeted therapy.

FIGURE 1 .
FIGURE 1. HFCWO device and imaging planes: (a) Typical modern HFCWO device with multiple actuators that can be activated individually for smart therapy, and (b) Anatomical imaging planes.

FIGURE 2 .
FIGURE 2. Model of human respiratory airway system: (a) airway tree of bifurcating segments, (b) transmission line model of the segment and its equivalent circuit with lumped parameters.

FIGURE 3 .
FIGURE 3. Model of respiratory airways by the equivalent circuit with lumped admittance parameters: (a) Node and branch indices, with encircled numbers representing the n number nodes in the branch order, and underlined numbers denoting the b branch order, and (b) standard branch.

FIGURE 4 .
FIGURE 4. The vertical and horizontal lines separate the airway geometry with the multiple sensing areas, and the known position of the simulated acoustic sensor array design is denoted with circles.

FIGURE 7 .
FIGURE 7. Quantitative lung function assessment through the mean image intensity and the thickening factor of AWT.

FIGURE 8 .
FIGURE 8. Acoustic imaging of obstructed airway with AWT increased by about a factor of 1.7.(a) , (c) Acoustic image produced from lung sound signal.(b), (d) Model acoustic image produced from airway pressure signal.

FIGURE 10 .
FIGURE 10.The relation between sensing sensitivity and the minimal nidus length that can be observed with a predetermined sensor number.

FIGURE 11 .
FIGURE 11.The relationship between sensor number and minimal nidus that can be observed on the right posterior of the chest wall.(a) The theoretical impact of sensor number required to identify the nidus length, and (b) The typical sensor numbers in a practical acoustic imaging system.

TABLE 1 .
Literature reviews key points.

TABLE 2 .
Incidence matrix, branch admittance vector, and branch voltage source vector.

TABLE 3 .
Material parameters of the airway geometry.

TABLE 4 .
Summary of key information about airway remodeling and airway wall thickness.

TABLE 5 .
Summary of key advantages and limitations of previous work and this study.