Rohan
Nalumachu‡
,
Anna
Anandita‡
and
Dharitri
Rath
Department of Chemical Engineering, Indian Institute of Technology, Jammu, Jammu and Kashmir-181221, India. E-mail: dharitri.rath@iitjammu.ac.in
First published on 15th March 2023
Competitive immunoassays are important diagnostic assays for the detection of small molecules such as vitamins, minerals, or some hormones. Although these assays are traditionally used to quantify small molecules, they are not extensively integrated with the paper-based devices. Numerical prototyping of these assays would be of paramount importance as it can help prior design of the devices, and therefore can reduce the time and resources needed. In this work, we are the first ones to present a thorough analysis of the computational model of the paper-based competitive immunoassay. The governing physics along with the pertinent boundary conditions coupled with the reactions both at the test line and control lines were considered to model this system. Furthermore, the performance of the device was evaluated through a simpler scaling analysis. Three important non-dimensional parameters were identified as T/C ratio, Pe, and Da to design such paper-based devices, and a design framework was presented to the readers. This opens up avenues for researchers to gain prior knowledge on their device performances.
Paper-based devices can be used to design diagnostic assays including immunoassays or nucleic-acid based assays.8 The immunoassays employed include either a sandwich format or a competitive format. The sandwich immunoassay format employs two antibodies that bind to distinct locations on the antigen and is frequently employed in the development of lateral flow devices, which is also the design employed in the most acknowledged paper-based device i.e., “pregnancy kits”. On the other hand, in the case of competitive assays, an antigen molecule and a labelled antigen compete for a limited number of antibody binding sites.9 Competitive assays have been primarily used for the detection of small molecules9 such as the detection of vitamins, minerals, and other small molecules.10 The development of competitive immunoassays on paper-based microfluidic substrates, therefore, would be of particular interest to the development of rapid diagnostic kits for the detection of small molecules.11 However, there are only a few studies reported on competitive assay formats; specifically, a complete numerical prototyping is lacking. In the literature, a few studies include understanding the kinetics of these assays, without demonstrating the details of how these occur in a lateral flow assay format.9,11–15 For instance, Sotnikov et al. reported an analytical model that explains the reaction schemes for sandwich and competitive immunoassays without considering the transport properties.14 They developed a mathematical ICA model to characterize immunological complex development and studied the influence of some of the model parameters, using a general analytical solution. Bishop et al. pointed out the importance of understanding the influence of parameters such as the flow characteristics, surface contacts, receptor–ligand kinetic interactions, and the presence of non-specific binding components while designing competitive immunoassays; however, they did not propose how to model these complexities during device design.16 In contrast, numerical prototyping of sandwich immunoassays using advection–diffusion–reaction models has been elaborated in the literature more extensively.17 Fu et al. compared modelling and experimental assays of a simple two-dimensional (2D) finite element model of convection, diffusion, and binding within a microchannel.15 Furthermore, Gasperino et al. reviewed the mechanistic principles underpinning the operations of lateral flow sandwich immunoassays, thereby describing the pertinent parameters influencing the functioning of these paper-based devices.18 Recently, Sathishkumar et al. reported a similar transport–reaction model for prototyping the commercially available pregnancy kits detailing the sandwich assay format in a lateral flow device, wherein they considered a constant convective velocity across the flow domain.19 Furthermore, there are convection–diffusion–reaction (CDR) models described earlier to predict the signal intensity of sandwich-type lateral flow assays under various design parameters.20,21 However, they did not look into the effect of flow in porous media.
In the current study, we are the first ones to present a complete numerical package of a competitive immunoassay on paper-based lateral flow devices. We integrate the transport and reaction phenomena in these membranes to describe the rate of capture of target analytes in the membrane. In one of the recent articles by one of the authors, the coupling of transport and reactions was already demonstrated to be useful in device design.7 In this article, the model is exclusively designed for quantifying the integration of fluid flow along with that of the transport and reactions simulated in real-time using a Multiphysics numerical package. The flow in the porous paper membranes is taken into consideration.7 The coupled physics describes an integrated approach to designing such devices starting from the delivery of target molecules to the quantification of the signals obtained after they react with the immobilized probe molecules present on the sensor surfaces. Finally, a scaling analysis of the device is done with the use of pertinent non-dimensional parameters, which helps in having an objective design for developers. Based on the above framework, a table guiding potential researchers is presented which can be used to have a prior idea about designing diagnostic devices with increasing sensitivities.
Fig. 1 depicts the architecture of a competitive immunoassay designed in a lateral flow format. It consists of basic components similar to any lateral flow assay device such as a sample pad, conjugate pad, test membrane and wicking pad. The conjugate pad is functionalized with reporter molecules (P), which in this case, are the primary antibodies conjugated to report nanoparticles (gold NPs in most cases). The sample pad is where a user adds the sample having target analyte (A) molecules. The test line is pre-functionalized with competitive antigen molecules (C) and the control line with secondary antibodies (Q). The reactions which take place on the lateral flow device operated in a competitive format are:
Bulk reaction:
A + P ↔ AP | (1) |
C + P ↔ PC | (2) |
Q + P ↔ QP | (3) |
Q + AP ↔ QPA | (4) |
QP + A ↔ QPA | (5) |
In order to simulate the above system, it's important to know the governing physics which involves the flow because of capillary pressure, coupled to the transport of species along with the various reactions as shown in Fig. 1. The governing physics accounts for the capillary flow and transport of species, and calculates the kinetics of both bulk and surface reactions by coupling the flow with transport. Eqn (6) presents the species transport along with the reactions solved for each species.
![]() | (6) |
![]() | (7) |
![]() | (8) |
Parameters | Notation | Value [units] | Ref. |
---|---|---|---|
Forward rate constants for bulk and test line reaction | KBf1, KTf1 | 1.97 × 101 m3 mol−1 s−1 | 23, 25 |
Forward rate constants for control line reactions | KCf1, KCf2, KCf3 | 1.97 × 101 m3 mol−1 s−1 | 23, 25 |
Backward rate constants for bulk and test line reaction | KBd1, KTd1 | 1.08 × 10−4 s−1 | 23, 25 |
Backward rate constants for control line reactions | KCd1, KCd2, KCd3 | 1.08 × 10−4 s−1 | 23, 25 |
Diffusion coefficients of A, P, and AP, respectively | D A, DP, DAP | 1× 10−10 m2 s−1 | 26 |
Inlet or initial concentrations of P, C, and Q | P in, Cin, Qin | 5 × 10−5, 4 × 10−3, 4 × 10−3 mol m−3 | 14, 27 |
Inlet concentration range of cortisol | A in | 1 × 10−6 to 5 × 10−4 mol m−3 | 28, 29 |
Density of sites available for surface reaction | Λ s | (Cin + Qin) × 0.00125 mol m−2 | |
Site occupancy number of all surface compounds | Σ © (for surface compounds) | 1 | |
Porosity | εP | 0.7 | 30 |
van Genuchten parameters | α, n, l (sample and conjugate pads, a) | 4.3 m−1, 1.41, 0.5 | 1 |
van Genuchten parameters | α, n, l (working membrane, b) | 1 m−1, 2.33, 0.5 | 1 |
van Genuchten parameters | α, n, l (wicking pad, c) | 4.3 m−1, 1.41, 0.5 | 1 |
Hydraulic conductivity (Ks) | K s | a) 2.08 × 10−5 | 1 |
b) 5.67 × 10−6 | |||
c) 2.08 × 10−3 [m s−1] | |||
Pressure | p 0 | −0.0002 [Pa] | 1 |
Length of the sample pad | L s | 2 [cm] | |
Distance from the conjugate pad to the test line and the control line to the wicking pad | Sps | 0.5 [cm] | |
Length of test and control lines | L th | 1.1 [mm] | |
Width of the LFIA | T s | 0.125 [cm] | |
Pressure head | H p | −5.3 [m] | 1 |
Maximum saturation | θ s | 1 |
Our model has a few assumptions; however, they do not affect the output signal measured from the device, and hence have been ignored. Our model ignores the effect of evaporation on the flow and transport and the dry storage of the report particles, and it assumes a uniform distribution of analytes in the sample pad. The evaporation effects can be neglected as the sample addition is really fast. The dry storage of the reporter particles in the conjugate pad could not be simulated. Instead, we considered that the report particles were added along with the analytes, i.e., the analytes and the reporters are added simultaneously to the sample pad, which then flow together. The model also ignores the uniform distribution of analytes in the sample pad. This could be a fair assumption as the time taken for uniform saturation of the sample pad is minimal.
![]() | (9) |
![]() | (10) |
As it can be observed for Fig. 4, there is a small variation between the literature and computed results and this can be attributed to the different numerical techniques used. Hence, the exclusive coefficient of determination (R2) values for each of the reported analyte concentrations were obtained (Fig. 4). The R2 values (regression coefficients) serve as an indication that the COMSOL model closely predicts the behaviour obtained through simulation results in the referenced literature. The regression analysis is as follows:
The R2 value determines how close two sets of values are. If R2 = 1, the simulation results are similar to that in referenced literature. However, some variation may be attributed to the use of different numerical techniques by our group and the authors of the studies. Our COMSOL model could predict the behaviour obtained by their simulation results closely, with evident coefficient of determination (R2, regression coefficient) values of ∼99%, showcasing the similarity. Further, it can be observed that the nature of the T/C ratio varies for the different concentration ranges because of the hook effect, which solely depends on the reactions that take place on the test line of the sandwich immunoassay. A detailed explanation of the behaviour is elaborated in ESI† S4 for the benefit of interested readers. Further, we also show a plot depicting the pattern of the signal formed on the control line for a wide range of analyte concentrations in ESI† S5 which substantiates the importance of test line signals in explaining the hook effect.
![]() | ||
Fig. 5 T/C vs. number of elements used for meshing, showing the optimized mesh used for designing the assay. |
Theoretically, by increasing the concentration of the competitor molecules, the magnitude of the test line signal should increase, which leads to an increase in the T/C value. But, from Fig. 6A-inset, it's conspicuous that the T/C value is greater for a competitor concentration corresponding to 4 × 10−3 (mol m−3) than 4 × 10−2 (mol m−3), for every time step, which insinuates a greater PC concentration for Cin = 4 × 10−3 (mol m−3) as mentioned in Fig. 6A-inset. But, here, the concentration of PC varies in accordance with Le Chatelier's principle: as the concentration of PC increases on the test line, the occurrence of backward reaction increases and results in a lower signal on the test line, indeed lowering the T/C value. And, for Cin = 4 × 10−4 (mol m−3), as the competitor concentration is relatively low, the PC formed will be lower and so will be the T/C value, as compared to the T/C values corresponding to both Cin = 4 × 10−2 (mol m−3) and 4 × 10−3 (mol m−3). The trend in Fig. 6A validates the above statements and implies Cin = 4 × 10−3 (mol m−3) as the suitable competitor concentration for this particular competitive LFIA to show a maximum signal on the test line. From Fig. 6B-inset, we can observe that the control line concentration corresponding to the initial concentration of secondary antibodies Qin = 4 × 10−3 (mol m−3) is greater than those corresponding to Qin = 4 × 10−4 (mol m−3) and even Qin = 4 × 10−2 (mol m−3), for all time steps. This means that a similar trend to that observed in the change of Cin regarding the test line concentration is observed here corresponding to the control line concentration. The trend in Fig. 6B upholds the above statements and implies Qin = 4 × 10−3 (mol m−3) as the suitable competitor concentration for this particular competitive LFIA to show a maximum signal on the control line.
The signal at the test line corresponds to the attachment of reporter molecules with the competitor molecules. The colorimetric signal corresponding to the formation of this complex PC is inversely related to the analyte A concentration, as it is the competitor molecule that gives a signal. The transport–reaction model along with the boundary conditions was solved to obtain the signals at the test and control lines. Fig. 7 illustrates the T/C ratio for various analyte concentrations, obtained with our proposed LFIA device used in this study. Fig. 7A shows a nonlinear behaviour of the signal obtained for the concentration ranges (1 × 10−6 to 5 × 10−4 mol m−3) of cortisol considering its presence in clinical samples corresponding to a pathological condition.28,29 Please note that the signal generated would have a contribution from both the test line and control line, wherein the signal generated on the test line would be inversely proportional to the analyte concentration as we are measuring the signal for a competitive assay in this case. The test line and control line signals were obtained separately from COMSOL simulations and plotted for our references in ESI† S7. One important observation for the control line signal graph was that the concentration of QPA + QP on the control line is not significantly affected by varying the concentration of analyte (ESI† S7), since either AP or P is abundant at the control line.3,7,21,32 In accordance with the above design, as observed in Fig. 7A, the signal corresponding to the lowest concentration (1 × 10−6 mol m−3) has the highest value, and vice versa. The signal intensity reached a maximum value at 250 s for the chosen set of conditions for cortisol detection, after which it decreases to reach the saturation value.
![]() | ||
Fig. 7 A. T/C vs. analyte concentration for different time steps. B. T/C vs. time for different analyte concentrations (1 × 10−6 to 5 × 10−4 mol m−3). |
In addition to the above observation, it was interesting to look at the kinetics of the signal intensity at different concentration values, which is presented in Fig. 7B. As shown in the figure, a non-monotonic behaviour is observed in the variation of T/C with time wherein the signal first increases and then decreases, for all the analyte concentrations considered in this work. The signals at the test line come from the formation of PC whereas they come from the formation of both QP and QPA at the control line. Now, let us categorize our concentrations into either low (black, red, blue and green lines in Fig. 7B) or high analyte concentrations (purple and yellow lines in Fig. 7B), as we can see slightly different trends for both these ranges. Now, low concentrations of A in the sample implies a higher concentration of unbound P available for reaction on the test line. As a result, the concentration of PC formed at the test line increases at a faster rate in the initial times till it reaches saturation. Meanwhile, the rate of formation of QPA at the control line occurs at a lower rate because of the following reasons: (a) the concentration of A is low, (b) the availability of unbound P from at the control line is lower as many of them have reacted at the test line itself, and (c) the formation of QP occurs at a lower rate because of a reduction in the concentration of Q due to the formation of QPA. Furthermore, it can be observed in the figure that at higher time (t ≥ 250 s in this case), after the signal reaches saturation, the availability of P on the test line with time decreases at the test line, and increases at the control line; hence the T/C decreases till the control line reaches equilibrium. In summary, though the concentrations of PC and QPA + QP increase, PC forms and saturates more rapidly than any of the control line reactions.
Also, the rate of reduction in the T/C ratio is steeper because the concentration of QPA increases and saturates at a lower rate. Similarly, for a higher analyte concentration where the concentration of unbound P drops, the quantity and rate of formation of PC on the test line is almost negligible, and the rate of formation of QPA increases, and so does the rate of saturation, on the control line. Hence, the T/C profiles for higher analyte concentrations are almost level, as shown in Fig. 7B.
Pe indicates the relative importance of convection and diffusion. In this study, we calculated the spatiotemporal variation of the convective velocity using the Richards equation. It is feasible to do so as the material properties such as the porosity, pore-size distribution, and hydraulic conductivity are known for different parts of the device. On the other hand, diffusion is a molecular phenomenon which is dependent on the physicochemical properties of the antigen and the solvent, including the radius of the antigen particle, viscosity, medium temperature, and the material from which the device is manufactured. Pe was thus calculated by substituting the average value of velocity throughout the timeline when the sample resides in the working membrane (t = 100–300 seconds), the characteristic length, and the diffusion coefficient in eqn (9).
Fig. 8 presents the non-dimensional signal intensity (T/C ratio) with the variation of Pe, for an inlet analyte concentration (Ain) of 10−6 mol m−3. Theoretically, the concentration signals corresponding to the test and control lines should decrease with an increase in Pe, as an increasing Pe corresponds to either an increase in convective velocity, or a decrease in species diffusivity. This may lead to the fact that the reactants do not get sufficient time to react before they flow past the reaction zone, which is the test line or the control line. Consequently, we will get a lower capture of the analytes, lowering the signal intensity. Fig. 8A–C exactly depict the same observation when we plot the corresponding signal concentrations at the test line, and control line for Pe varying from 5 × 10−3 to 5.5 × 107. These plots were obtained for different Pe regimes, corresponding to the dominance of diffusion or convection. Firstly, considering a diffusion dominant regime (Pe < 300), as Fig. 8A shows, the decrease in control line concentration is greater in comparison to that of concentration on the test line. This could be attributed to the fact that the signal at the control line requires the diffusion of all three molecules where A, P, and AP have to concentrate on the control line to form a signal, as opposed to the test line which needs the availability of P alone. Hence, we plot the T/C ratios for this Pe range, and as shown in Fig. 8D, we observe a monotonic increase in the signal. As shown in Fig. 8B, the slope of concentration on the test line is steeper in comparison to that of concentration on the control line. As the sample velocities are low, the concentrations on the test and control lines depend on the reaction kinetics at the corresponding locations, i.e., the rate of PC and PQ + QPA formation on the test and control lines, respectively. Thus, the varied rates of reaction are portrayed in Fig. 8B. Next, we observe a non-monotonic variation of T/C with an increase in Pe (300 < Pe < 5.5 × 107) as shown in Fig. 8E. Here, it can be observed that there is a sudden decrease in the T/C signal corresponding to Pe = 275. This is attributed to the fact that at this value of Pe, there is a transition from a diffusion dominant to a convection dominant regime. As the sample velocities are low (550 ≤ Pe < 5.5 × 104), the concentrations on the test and control lines depend on the reaction kinetics at the corresponding locations, i.e., the rate of PC and PQ + QPA formation on the test and control lines, respectively. Thus, the varied rates of reaction are portrayed in Fig. 8B. Therefore, Fig. 8E demonstrates that the range of 550 ≤ Pe < 5.5 × 104 led to a decrease in T/C. Considering Pe ≥ 5.5 × 104, approximately the same slopes were obtained for the test and control line signals, as shown in Fig. 8C, which is explained by the non-dependence on reaction rates. The observed effect of Pe ≥ 5.5 × 104 is an increase in T/C, which correlates positively with an increase in Pe. This relationship is depicted in Fig. 8E as a linear slope with an angle of approximately 45 degrees. Hence, showing a non-monotonic behaviour of T/C with varying Pe (operated convective dominant regime), and raising a rather inquisitive nature for variation between rate of convection and reaction. Therefore, the Damköhler number (Da) is used to explicitly portray the relevance of kinetics to the signal (T/C).
![]() | ||
Fig. 8 Signal intensities on the test line (P), and control line (QP + QPA) vs. Pe for A. Pe < 10, B. 10 ≤ Pe < 1000, and C. Pe ≥ 1000. T/C ratio vs. Pe for D. Pe < 10 and E. Pe ≥ 10. While computing, the LFIA analyte concentration is taken as 10−6 mol m−3, and the saturated values of concentration on the test and control lines and the T/C ratio are plotted. The calculations specific to this device detecting cortisol are given in ESI† S1. |
The Damköhler number (Da) is an essential consideration for any system wherein contributions from both transport and reactions are important. Da compares the rates of flow and reaction kinetics and Da is calculated using eqn (10). Fig. 9 depicts the variation of T/C with respect to Da. For higher Da (Da = 600), the kinetics dominate over the transport. Therefore, the signal becomes transport limited. For values of Da < 600, T/C is dependent on the convection. Hence in this regime, as Da increases, say from 15 to 300, the convective velocity decreases, leaving more time for the reaction to happen. As a result, as presented in Fig. 9A illustrates an increase in signal on the test line. However, the faster reaction rate at the test line leads to a higher rate of consumption of reporter particles P, leaving less of them that reach the control line. This explains the decrease in the signal at the control line shown in Fig. 9B. Therefore, the combined effect is reflected in Fig. 9C for an increase in Da (for Da < 500). Now, we consider the regime for which Da > 600. Theoretically, an increase in Da implies an increase in reaction kinetics, which thereby should increase the test line concentration and indeed the T/C value. However, this is not the trend we obtain in Fig. 9A–C. The anomaly could be explained using the Le Chatelier's principle. An increase in the forward kinetics of reaction on the test line (bubble 2 in Fig. 1) results in an increase in the concentration of PC, thereby increasing the occurrence of backward reaction and lowering the signal on the test line, and lowering the value of T/C. However, Fig. 9B shows that the signal at the control line is not affected much with the reversibility as the consumption of P depends on three simultaneous reactions. Further, Fig. 9C shows the combined effect from Fig. 9A and B, wherein the T/C ratio happens to decrease when Da is increased (for Da ≥ 600).
![]() | ||
Fig. 9 Signal intensities on the A. test line (P) and B. control line (QP + QPA) vs. Da. C. T/C vs. Da for 0.01 ≤ Da ≤ 100. |
The pertinent input and output parameters that would guide researchers to design any such lateral flow based competitive immunoassay platform are hence summarised in Fig. 10. Briefly, the design of such types of immunoassays includes first selecting the input parameters starting with choosing the right kinetics pertinent to the reactions, then obtaining the material and transport properties to solve the flow and transport equations, and finally, the use of scaling analysis to obtain the regime of operations and a decision on the final output parameters.
Footnotes |
† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d2sd00211f |
‡ Equal contribution by the authors. |
This journal is © The Royal Society of Chemistry 2023 |