Expediting Disulﬁram Assays through a Systematic Analytical Quality by Design Approach

: An Analytical Quality by Design (AQbD) approach is presented, aiming at the development and validation of an HPLC method for the quantiﬁcation of disulﬁram and copper diethyldithiocarbamate in lipid nanoparticles. Following the deﬁnition of the analytical target proﬁle (ATP), encompassing the critical analytical attributes (CAA), a two-level risk assessment strategy (Ishikawa diagram—failure mode and effect analysis (FMEA)) was employed to identify the critical method parameters (CMPs) with an extensive impact on method performance. The behavior of the CMPs (ﬂow rate and mobile phase composition) was further characterized by experimental design, resorting to a face-centered central composite design (FcCCD). Statistical modeling, response surface analysis, and Monte Carlo simulations led to the deﬁnition of the Method Operable Design Region (MODR), associated with a negligible risk of failing the predeﬁned CAA speciﬁcations. The optimal method was validated according to international regulatory recommendations. Apart from guaranteeing linearity, accuracy, precision, speciﬁcity, robustness, and stability, these conditions were found to be suitable for analysis using a different HPLC column and equipment. In a nutshell, the development and optimization strategies, under the comprehensive framework of AQbD, provided an effective, simple, rapid, reliable, and ﬂexible method for routine analysis of the compounds in research or industrial environments.


Introduction
Disulfiram (DSF) is a dithiocarbamate derivative with clinical application to treat alcohol addiction. It blocks ethanol metabolism by inhibiting hepatic aldehyde dehydrogenases 1 and 2, thus increasing acetaldehyde blood levels. Consequently, alcohol consumers experience nausea, sweating, hypotension, respiratory difficulties, and other alcoholic intoxication symptoms, thereby rejecting additional beverages [1]. In parallel, DSF shows promising in vitro/in vivo results against several types of neoplastic diseases, as it may act over 19 different targets/pathways to reduce cancer cell viability [2]. DSF is also a chelating compound, producing copper (II) diethyldithiocarbamate (Cu(DDC) 2 ) in the presence of endogenous or exogenous copper (II) ions. Alternatively, Cu(DDC) 2 can be synthesized by chelating copper with sodium diethyldithiocarbamate (Figure 1) [3,4]. Interestingly, this However, there is a lack of clinical evidence that supports the efficacy of DSF as an anticancer drug due to the extensive degradation and/or metabolism (reduced half-time) of the drug. The reduced aqueous solubility of DSF (0.2 g/L) also impairs its intravenous administration. As for Cu(DDC)2, its intrinsic solubility is even lower, estimated as 0.700 µg/L [5]. Their inclusion in nanosystems for cancer drug delivery is thus highly anticipated.
The increasing interest in these compounds prompts the need for developing appropriate analytical methods, with chromatography remaining the standard approach for quantification in the pharmaceutical industry. Curiously, neither DSF nor Cu(DDC)2 possess ionizable groups that can modify their partition between the stationary and the mobile phases during chromatography. There are also no reports in the literature for their simultaneous quantification by reversed-phase high-performance liquid chromatography (RP-HPLC). An exception is presented by Irth et al. [6], dated from 1988, who used postcolumn derivatization for spectrophotometric analysis. Furthermore, reported chromatographic methods for DSF or Cu(DDC)2 focus on the quantification of the compounds in nanotechnological platforms, without indication of Analytical Quality by Design (AQbD) development or formal validation.
Similarly to Quality by Design (QbD), one can say that AQbD major goal is an indepth understanding of the method. Additionally, the method is well understood when all sources of critical variability are identified and characterized, with the method operable design region (MODR) accurately and successfully predicting the critical analytical attributes (CAAs) [7]. This translation to analytical procedures gained increasing attention over the last decade. Contrary to QbD, it is still scarcely evidenced in the literature, with only 86 references in the period of 2010 to 2018 [8].
AQbD includes the identification of the Analytical Target Profile (ATP) and CAAs, as well as the critical method parameters (CMPs) that significantly affect method performance [9]. The identification of these parameters is often conducted through risk assessment, either using Ishikawa diagrams, risk matrices, failure mode and effect analysis (FMEA), or failure mode, effects, and criticality analysis (FMECA), although other tools can be used [10]. In every case, there is a strong contribution from previous knowledge and experience. CMPs are then systematically explored by Design of Experiments (DoE), an AQbD statistical tool for mathematical modeling that characterizes the factor-response (CMPs-CAAs) relationships through simultaneous variations of the CMP conditions. When several CMPs are identified, one should consider two-level designs, such as 2 k full factorial designs (FFD), fractional factorial designs (FrFD), and Placket-Burman designs (PBD). Moreover, also known as screening designs, these have proved their usefulness in confirming the criticality of the CMPs, often leading to their disregard for further studies. The three-level designs (response surface methodology), including 3 k FFDs, Box-Behnken designs (BBD), and central composite designs (CCD), are However, there is a lack of clinical evidence that supports the efficacy of DSF as an anticancer drug due to the extensive degradation and/or metabolism (reduced half-time) of the drug. The reduced aqueous solubility of DSF (0.2 g/L) also impairs its intravenous administration. As for Cu(DDC) 2 , its intrinsic solubility is even lower, estimated as 0.700 µg/L [5]. Their inclusion in nanosystems for cancer drug delivery is thus highly anticipated.
The increasing interest in these compounds prompts the need for developing appropriate analytical methods, with chromatography remaining the standard approach for quantification in the pharmaceutical industry. Curiously, neither DSF nor Cu(DDC) 2 possess ionizable groups that can modify their partition between the stationary and the mobile phases during chromatography. There are also no reports in the literature for their simultaneous quantification by reversed-phase high-performance liquid chromatography (RP-HPLC). An exception is presented by Irth et al. [6], dated from 1988, who used post-column derivatization for spectrophotometric analysis. Furthermore, reported chromatographic methods for DSF or Cu(DDC) 2 focus on the quantification of the compounds in nanotechnological platforms, without indication of Analytical Quality by Design (AQbD) development or formal validation.
Similarly to Quality by Design (QbD), one can say that AQbD major goal is an in-depth understanding of the method. Additionally, the method is well understood when all sources of critical variability are identified and characterized, with the method operable design region (MODR) accurately and successfully predicting the critical analytical attributes (CAAs) [7]. This translation to analytical procedures gained increasing attention over the last decade. Contrary to QbD, it is still scarcely evidenced in the literature, with only 86 references in the period of 2010 to 2018 [8].
AQbD includes the identification of the Analytical Target Profile (ATP) and CAAs, as well as the critical method parameters (CMPs) that significantly affect method performance [9]. The identification of these parameters is often conducted through risk assessment, either using Ishikawa diagrams, risk matrices, failure mode and effect analysis (FMEA), or failure mode, effects, and criticality analysis (FMECA), although other tools can be used [10]. In every case, there is a strong contribution from previous knowledge and experience. CMPs are then systematically explored by Design of Experiments (DoE), an AQbD statistical tool for mathematical modeling that characterizes the factor-response (CMPs-CAAs) relationships through simultaneous variations of the CMP conditions. When several CMPs are identified, one should consider two-level designs, such as 2 k full factorial designs (FFD), fractional factorial designs (FrFD), and Placket-Burman designs (PBD). Moreover, also known as screening designs, these have proved their usefulness in confirming the criticality of the CMPs, often leading to their disregard for further studies. The three-level designs (response surface methodology), including 3 k FFDs, Box-Behnken designs (BBD), and central composite designs (CCD), are used for method optimization. Other two-and three-level designs can be used under the AQbD umbrella for method development and optimization, including Taguchi, I-optimal, D-optimal, and Doehlert designs [11].
Ultimately, AQbD defines a multidimensional region of robustness within the knowledge space, the MODR, in which the chromatographic conditions (in terms of CMPs) satisfy the criteria defined in the ATP [12,13]. This region can be constructed by overlaying the independent 2D contour plots that model CAA parameters in function of the variation of the CMPs. Nonetheless, Bayesian and stochastic approaches, supported by the determination of the risk of failing the specifications, have demonstrated superiority [8]. The MODR is also subject to validation, often by Monte Carlo simulations and/or experimental validation, in order to guarantee its authenticity and robustness. The definition of the MODR encompasses various advantages. Apart from increasing method knowledge, it indicates a range of values in which the method provides quality results, thus being suitable for analysis. In fact, the method parameters can be varied anywhere within the MODR, in opposition to a fixed condition, without the need for revalidation, but only for an adjustment [8,14]. Note that an appropriate control strategy must be implemented to control all sources of variability that might occur, always aiming at the continual improvement of the method [14]. An overview of the AQbD workflow is presented in Figure 2.
Chemosensors 2021, 9, x FOR PEER REVIEW 3 of 20 used for method optimization. Other two-and three-level designs can be used under the AQbD umbrella for method development and optimization, including Taguchi, I-optimal, D-optimal, and Doehlert designs [11]. Ultimately, AQbD defines a multidimensional region of robustness within the knowledge space, the MODR, in which the chromatographic conditions (in terms of CMPs) satisfy the criteria defined in the ATP [12,13]. This region can be constructed by overlaying the independent 2D contour plots that model CAA parameters in function of the variation of the CMPs. Nonetheless, Bayesian and stochastic approaches, supported by the determination of the risk of failing the specifications, have demonstrated superiority [8]. The MODR is also subject to validation, often by Monte Carlo simulations and/or experimental validation, in order to guarantee its authenticity and robustness. The definition of the MODR encompasses various advantages. Apart from increasing method knowledge, it indicates a range of values in which the method provides quality results, thus being suitable for analysis. In fact, the method parameters can be varied anywhere within the MODR, in opposition to a fixed condition, without the need for revalidation, but only for an adjustment [8,14]. Note that an appropriate control strategy must be implemented to control all sources of variability that might occur, always aiming at the continual improvement of the method [14]. An overview of the AQbD workflow is presented in Figure 2. Here, an AQbD strategy is set for the quantification of DSF and its anticancer metabolite, Cu(DDC)2, in nanostructured lipid carriers (NLCs), the second generation of solid lipid nanoparticles.

Copper (II) Diethyldithiocarbamate Synthesis
The synthesis of Cu(DDC) 2 was achieved through the chelation of copper ions by sodium diethyldithiocarbamate in an aqueous environment, as described by Liu et al., with modifications [15]. Briefly, copper (II) sulfate (1 mmol, 0.253 g) and sodium diethyldithiocarbamate (2 mmol, 0.456 g) were mixed in a 1:2 molar ratio, under stirring, at room temperature, for 1 h. Following the spontaneous reaction, the solid brown precipitate was filtered under reduced pressure, washed with ultrapurified water, recrystallized with ethanol to eliminate salts and unreacted products, and dried under vacuum (91% yield). Structural analysis was confirmed by 1 H and 13 C NMR (Supplementary Material, Figures S1 and S2).

Chromatographic Instrument and Conditions
The chromatographic analysis was carried out in a Shimadzu LC-2010HT system, equipped with a quaternary pump (LC-20AD), an autosampler (SIL-20AHT), an oven (CTO-10AS), and a photodiode array detector (SPD-M2OA). The separation was achieved utilizing a reversed-phase ACE ® 5 C 18 column, with 5 µm of particle size, 4.6 mm of internal diameter, and 150 mm of length (Advanced Chromatography Technologies Ltd., UK), supported by a SecurityGuard cartridge under the following conditions: isocratic mode, at 40 • C, with a mobile phase composed of ACN:H 2 O in various ratios and different flow rates. Using an injection volume of 10 µL, DSF and Cu(DDC) 2 were determined at their maximum absorption wavelength, 217 nm and 270 nm, respectively.

Calibration Standards and Quality Control Solutions
DSF (5 mg in 5 mL) and Cu(DDC) 2 (2.5 mg in 5 mL) stock solutions were prepared independently in ACN and sonicated to guarantee a complete dissolution. Note that methanol, tetrahydrofuran, dimethyl sulfoxide, and dichloromethane led to an accelerated degradation of the compounds or showed insufficient solvent abilities. Two working standards containing both compounds at 100 and 10 µg/mL were prepared by dilution of the stock solutions. Calibration standards (0.1, 0.25, 0.5, 1, 5, 10, 25, 50, 75, and 100 µg/mL) were prepared by sequential dilutions of the working solutions with mobile phase and filtered through 0.2 µm PTFE filters prior to the injection. Calibration curves were constructed by linear regression of the peak area against the nominal concentration after weighting correction. Quality controls were prepared similarly to calibration standards, at intermediate concentrations (0.2, 20, and 90 µg/mL).

Analytical Target Profile Definition
The establishment of an analytical QbD approach assumes the establishment of an analytical target profile (ATP), a prospective summary of the quality characteristics of the method, similarly to the quality target product profile (QTPP) in pharmaceutical development. Table 1 indicates various ATP elements taken into consideration for the development and optimization of the chromatographic method for the quantification of DSF and its anticancer metabolite, Cu(DDC) 2 . Considering their relevance in liquid chromatography, the number of theoretical plates (N), retention time (R t ), tailing factor (T f ), critical peak resolution (Res), and capacity factor (k') of the analytes were selected as critical analytical attributes (CAAs), and should be consistent with formal or commonly acceptable validation criteria.

Risk Assessment
Risk assessment approaches are useful to identify potential analytical parameters, such as method parameters, equipment, and measurements that may have an impact on method development and performance. The simplest approach is the conceptualization of Ishikawa (also known as fishbone and cause-and-effect) diagrams that outline the relationship between possible critical method parameters and the critical analytical attributes established in the ATP. Based on the analytical parameters identified in the Ishikawa diagram, a failure mode and effect analysis (FMEA) was used to evaluate and select the critical ones to be further studied. Risk priority number (RPN) was determined based on, (1) where an RPN above 100 was considered of high risk and failure mode and therefore submitted to a systematic evaluation using a response surface method optimization [16]. The criteria considered for the determination of the RPN are displayed in Table 2.

ATP Element Target Justification
Analyte DSF, Cu(DDC) 2 Development of an analytical method for the estimation of the analytes in solution or co-encapsulated in nanosystems for routine and stability analysis.

Sample Liquid
Analytes must be suitably dissolved, with complete solvent miscibility. Analyticaltechnique

RP-HPLC
Highly lipophilic analytes are better separated and eluted in RP-HPLC methods due to the increased retention with the non-polar C 18 stationary phase. Apart from the increased resolution, RP-HPLC provides a fast analysis with small sample volume and organic solvent consumption.

Instrument
HPLC with a quaternary pump, autosampler, and PDA detector A quaternary pump allows an improved mixing of the mobile phase solvents and higher resolution, whereas the PDA detector allows the detection of several compounds at their λ max , thus maximizing sensitivity (DSF, λ max = 217 nm, Cu(DDC) 2 , λ max = 270 nm).

Method
Specific, accurate, precise, linear, reproducible, robust, cost-effective, and simple The method should comply with formal validation criteria, presenting a short run time and reduced use of organic solvents.

Application Assay
Method should be able to determine DSF and Cu(DDC) 2 in solution or co-encapsulated in nanosystems for routine and stability analysis.

CAAs
Number of theoretical plates, retention time, tailing factor, critical peak resolution, capacity factor Chromatographic parameters that allow a robust and reliable determination of the analytes. Should meet their formal and commonly accepted quality criteria.

Method Development and Optimization
The chromatographic method was developed and optimized using a face-centered central composite design (FcCCD), an effective second-order model that combines a twolevel factorial design (2 n ) or fractional factorial design (2 n−k ), one central point and 2n outer points, denoted as axial or star points (along with the coordinate axes, at ± α). In the case of an FcCCD, with two factors, both factorial and star points present the same negative and positive distance from the central point, i.e., the star points are at the center of each face of the experimental domain (α = 1), at (±α, 0) and (0, ±α), thus being identical to a 3 2 FD with a rectangular experimental domain ( Figure 3) [17]. central composite design (FcCCD), an effective second-order model that combines a twolevel factorial design (2 n ) or fractional factorial design (2 n−k ), one central point and 2n outer points, denoted as axial or star points (along with the coordinate axes, at ± α). In the case of an FcCCD, with two factors, both factorial and star points present the same negative and positive distance from the central point, i.e., the star points are at the center of each face of the experimental domain (α = 1), at (±α, 0) and (0, ±α), thus being identical to a 3 2 FD with a rectangular experimental domain ( Figure 3) [17]. The proposed models for each CAA are explained by the quadratic polynomial function: where Y is the evaluated response associated with each factor level combination, β0 is the response in the absence of effects, β1 and β2 are the linear coefficients of the factors X1 and X2, respectively, β12 is the interaction coefficient between X1 and X2, and β11 and β22 are the quadratic terms that allow the prediction of the curvature of the model. The statistical significance of the parameters in the regression model was evaluated using Student's t-test (95% confidence level, α = 0.05) and ANOVA (p-value < 0.05). Non-significant terms were removed by backward selection, in which higher p-value terms are sequentially removed (p-value > 0.05 threshold). FcCCD models were built and analyzed in JMP ® Pro 15.0.0. The proposed models for each CAA are explained by the quadratic polynomial function:

MODR Design and Validation
where Y is the evaluated response associated with each factor level combination, β 0 is the response in the absence of effects, β 1 and β 2 are the linear coefficients of the factors X 1 and X 2 , respectively, β 12 is the interaction coefficient between X 1 and X 2 , and β 11 and β 22 are the quadratic terms that allow the prediction of the curvature of the model. The statistical significance of the parameters in the regression model was evaluated using Student's t-test (95% confidence level, α = 0.05) and ANOVA (p-value < 0.05). Non-significant terms were removed by backward selection, in which higher p-value terms are sequentially removed (p-value > 0.05 threshold). FcCCD models were built and analyzed in JMP ® Pro 15.0.0.

MODR Design and Validation
The MODR was designed using Monte Carlo simulations, considering the specifications of each CAA. This operable region provides several conditions that fulfill the chromatographic requirements set in the ATP, expressed as the probability of failure (%), assuming a normal distribution of the method parameters associated with a certain random standard deviation [18]. Random variation for each CAA was also considered in a total of 10,000 experimental runs across the range of parameters. In addition, experimental validation, using 4 additional chromatographic conditions was also performed.

Method Validation
The validation of the optimized chromatographic method was performed according to ICH and FDA regulations, taking into account the system suitability, limits of detection and quantification, linearity, accuracy, precision, robustness, ruggedness, specificity, and stability [19][20][21][22].

System Suitability
System suitability parameters of the method were determined by injecting six times the same quality control solutions and compared to the specifications and recommendations of the FDA [21]. Accordingly, theoretical plate number (N) > 2000, tailing factor (T f ) < 2.0, critical peak resolution (Res) > 2.0, relative standard deviation (RSD) of peak area and of retention time (R t ) < 2.0% and capacity factor (k') > 2.0.

Detection and Quantification Limits
The limits of detection (LOD) and quantification (LOQ) were determined with a calibration curve obtained from specific calibration standards (0.1, 0.25, 0.5, 1, and 5 µg/mL). Subsequently, LOD and LOQ were calculated according to LOD = 3.3 σ/S and LOQ = 10 σ/S, where σ is the standard deviation of the response and S the calibration curve slope [19].

Linearity
The linearity of the method was evaluated using six calibration curves containing DSF and Cu(DDC) 2 , simultaneously ranging from 0.1 to 100 µg/mL. Linearity was determined through the plotting of peak area as a function of the nominal concentration of the standards, using a least squares regression method. Since homoscedasticity (condition of equal variance throughout the points) is not frequently observed, in particular, using wide xvalue ranges, a 1/x 2 weighting factor was employed (bias < 15%) [23]. According to regulatory guidelines, a determination coefficient, R 2 > 0.999 was considered [21].

Accuracy and Precision
Accuracy and precision were evaluated intra and interdaily for three consecutive days. Six independent samples of each quality control solution were prepared and analyzed. To demonstrate accuracy, i.e., the agreement between the nominal concentration (true) and the measured value, a limit of bias < 15% was established. As for intraday (repeatability) and interday (intermediate) precision, i.e., the degree of scatter among several determinations of the same homogeneous sample, a limit of RSD < 15% was imposed [22].

Robustness and Ruggedness
Robustness, as per the ICH definition, is a "measure of (the analytical procedure) capacity to remain unaffected by small, but deliberate variations in method parameters and provides an indication of its reliability during normal usage" [19]. As per that definition and validation of the MODR, the compliance to this parameter is empirically met. According to the USP, ruggedness is seen as a measure of the reproducibility of the results under various conditions met when changing the laboratory and the analyst [24]. Therefore, six quality control solutions were analyzed, for accuracy and precision, in a different HPLC column (Kinetex ® EVO C 18 column, with 5 µm of particle size, 4.6 mm of internal diameter and 150 mm of length) and HPLC instrument (Shimadzu LC-10AD, quaternary pump LC10-AD, autosampler unit Sil-10ADVP, CTO-10AVP oven and a CBM-20A UV detector). System suitability, according to Section 2.6.1 was also evaluated.
2.6.6. Specificity Specificity, termed as the ability to accurately measure the analytes in the presence of potential interferents, was visually analyzed in terms of peak R t [19]. A representative chromatogram of a quality standard solution and of the supernatant of lipid nanoparticles was compared to the supernatant obtained from a blank formulation [25]. The main goal of specificity is to ensure that the integrity of each analyte is not compromised by any formulation excipient (including the combination of lipids and surfactants).

Stability
In order to evaluate stability, quality control samples were analyzed for short-term stability, following 24 h at room temperature (approx. 25 • C, in the autosampler) and 72 h at 4 • C. Furthermore, stock solutions were also evaluated for three freeze-thaw cycles and intermediate-term stability, at −20 • C, for 7 days, by preparing and evaluating the respective quality control solutions [21,22].

NLC Production
NLCs were prepared by high shear homogenization followed by high-pressure homogenization, as described elsewhere [26]. Briefly, 2% (w/w) DSF was dissolved in the molten lipid phase, containing the solid and the liquid lipids (mono-, di-and triglyceride esters of fatty acids (C 10 to C 18 ), oleic acid and propylene glycol mono-and diesters of caprylic acid) as well as the dissolved salt of copper. Note that Cu(DDC) 2 was formed in situ, during the production method, through the reaction of DSF and copper ions. Subsequently, the mixture was added to 30 mL of a hot aqueous surfactant solution and emulsified for 5 min at 19,000 rpm (Ultra-Turrax X 10/25, Ystral GmbH, Dottingen, Germany). The pre-emulsion was further processed in a high-pressure homogenizer for 10 min at 1500 bar (150 MPa) (EmulsiFlex-C3; Avestin, Inc., Ottawa, Canada). The resulting dispersion was then stored at 4 • C to form the NLCs. Colloidal properties, in particular, particle size and distribution, and zeta potential, were determined following a 100-fold dilution in ultrapurified water by dynamic and electrophoretic light scattering (Zetasizer Nano ZS, Malvern Instruments, Malvern, UK).

Entrapment Efficiency and Drug Loading Determination
The entrapment efficiency (EE) and drug loading (DL) of DSF and CU(DDC) 2 in the NLCs were calculated according to where W total drug is the amount of drugs encapsulated in the NLCs, W free drug is the amount of free drugs determined in the aqueous phase after ultrafiltration-centrifugation (Sartorius ® Vivaspin 500 filter unit, 100 kDa molecular weight cut-off), and W lipid is the weight of the lipid phase of the nanosystems [25]. For the total drug quantification, a specific volume of NLCs was suitably diluted in the mobile phase and heated at 60 • C, under sonication, for 15 min. The solution was further centrifuged for 5 min at 11,740× g, at 4 • C, and the supernatant was filtered and analyzed by HPLC. For the determination of free drugs, NLCs were diluted with cold ACN in order to dissolve potential drug crystals, submitted to ultrafiltration-centrifugation, and the collected aqueous phase analyzed by HPLC.

Risk Assessment
Narrowing down to a few most critical method parameters for experimental evaluation can be challenging, considering the multitude of factors that potentially affect the development of a chromatographic method. A systematic risk assessment analysis not only identifies CMPs but also prioritizes them according to their severity. Here, a two-stage study was conducted, contemplating an Ishikawa diagram followed by a failure mode effect analysis (FMEA). The Ishikawa diagram was designed based on prior knowledge, both experimental and theoretical. This diagram outlines several factors that may contribute to method performance, in particular, considering the analyst, environment, equipment, method and measurement, and data ( Figure 4). Nonetheless, it does not consider their risk or criticality. development of a chromatographic method. A systematic risk assessment analysis not only identifies CMPs but also prioritizes them according to their severity. Here, a twostage study was conducted, contemplating an Ishikawa diagram followed by a failure mode effect analysis (FMEA). The Ishikawa diagram was designed based on prior knowledge, both experimental and theoretical. This diagram outlines several factors that may contribute to method performance, in particular, considering the analyst, environment, equipment, method and measurement, and data ( Figure 4). Nonetheless, it does not consider their risk or criticality. In order to efficiently categorize and identify the critical method parameters to be systematically studied, an FMEA approach was conducted. FMEA focuses on minimizing potential failures prior to their occurrence and aids in developing a strategy to reduce the risk of a future failure by changing the process appropriately [27]. Several factors identified in the Ishikawa diagram were evaluated for their severity, occurrence, and detectability (Table 3). Accordingly, flow rate and mobile phase composition were the method parameters with the highest RPN score, determined using Equation (1). The remaining failure modes, with ranging severity, occurrence, and detectability, were considered to be non-critical and can be easily resolved by the corresponding strategy. For example, line purging-related risks can be tackled by guaranteeing purge prior to the analysis, whereas UV detection (with respect to λmax, signal/noise ratio, and sensitivity) can be defined a priori, in the ATP.  In order to efficiently categorize and identify the critical method parameters to be systematically studied, an FMEA approach was conducted. FMEA focuses on minimizing potential failures prior to their occurrence and aids in developing a strategy to reduce the risk of a future failure by changing the process appropriately [27]. Several factors identified in the Ishikawa diagram were evaluated for their severity, occurrence, and detectability (Table 3). Accordingly, flow rate and mobile phase composition were the method parameters with the highest RPN score, determined using Equation (1). The remaining failure modes, with ranging severity, occurrence, and detectability, were considered to be non-critical and can be easily resolved by the corresponding strategy. For example, line purging-related risks can be tackled by guaranteeing purge prior to the analysis, whereas UV detection (with respect to λ max , signal/noise ratio, and sensitivity) can be defined a priori, in the ATP.

Method Development and Optimization: Critical Analytical Attributes
Taking into consideration the RPN > 100 threshold, flow rate (X 1 ) and mobile phase ratio (in % of ACN) (X 2 ) were defined as the CMPs and, consequently, systematically evaluated following the design matrix of Table 4. Their impact on CAAs (N, T f , Res, R t , and k') was characterized. The chromatographic conditions were set after preliminary studies, where other mobile phase compositions (in particular, with methanol and various pH values) were evaluated at different flow rates. Table 4. Face-centered central composite design (FcCCD) matrix for the evaluation of flow rate (X 1 ) and mobile phase composition (X 2 ). For each condition, two calibration standard solutions were considered (1 and 100 µg/mL) in order to inspect the impact of impurities on the critical analytical attributes (CAAs). According to the FcCCD, a total of 18 experimental runs were performed (nine chromatographic conditions, with two calibration standards, 1 and 100 µg/mL), at three levels: for X 1 , 0.8, 1.0 and 1.2 mL/min; for X 2 , 55, 70 and 85% of ACN. The influence of each factor and their combination in the CAAs was determined using the polynomial coefficients of Equation (2). Higher coefficient magnitudes indicate stronger main effects on the corresponding CAA. Additionally, the signal of the coefficient also has an important significance. When the coefficient sign is negative, an increase in the factor decreases the response. On the contrary, if the sign is positive, an increase in the factor increases the response. The interaction coefficient indicates how the variation of one factor modulates the effect of the other factor. Note that the main effect may be impacted by the interaction terms. Lastly, the quadratic terms, whenever statistically significant, provide curvature to the models [25]. Note that non-significant coefficients were removed by backward selection, as mentioned in Section 2.5.3. Therefore, only statistically significant terms (p-value < 0.05) are presented in the subsequent sections. The ANOVA parameters for the characterization of the model fitting per CAA are presented in Supplementary Material, Table S1.

Theoretical Plates
Theoretical plates are a measure of column efficiency, indicating the number of peaks located per run-time unit on a chromatogram [28]. Regardless of the conditions evaluated under the experimental design, N was always above the minimum required by regulatory guidelines (N > 2000) and can be estimated using: N DSF = 7390 − 939 X 1 − 2637 X 2 + 951 X 22 2 , R 2 = 0.946, R 2 adj = 0.935 N Cu(DDC)2 = 10,768 − 985 X 1 − 3484 X 2 , R 2 = 0.903, R 2 adj = 0.890 Considering the coefficients of Equations (5) and (6), both factors present a negative impact on this response, although with a higher extent for the % ACN (X 2 ). In addition, the β 22 coefficient provides curvature to the DSF model. In fact, as graphically evidenced by Figure 5, an increase in flow rate (X 1 ) and/or the percentage of the organic phase (X 2 ) reduces the number of theoretical plates for both analytes. This trend is also verified for other compounds under different chromatographic conditions [25]. NDSF = 7 390 − 939 X1 − 2 637 X2 + 951 X22 2 , R 2 = 0.946, R 2 adj = 0.935 (5) NCu(DDC)2 = 10 768 − 985 X1 − 3 484 X2, R 2 = 0.903, R 2 adj = 0.890 (6) Considering the coefficients of Equations (5) and (6), both factors present a negative impact on this response, although with a higher extent for the % ACN (X2). In addition, the β22 coefficient provides curvature to the DSF model. In fact, as graphically evidenced by Figure 5, an increase in flow rate (X1) and/or the percentage of the organic phase (X2) reduces the number of theoretical plates for both analytes. This trend is also verified for other compounds under different chromatographic conditions [25].

Tailing Factor
Chromatographic peaks are ideally represented by Gaussian curves. However, in practice, they are usually asymmetrical due to the interaction of the analytes and the stationary phase. Thus, Tf is a measure of peak symmetry [28]. According to official guidelines, and to reduce peak integration variability, Tf < 2.0. For each analyte, Tf can be determined by Tf DSF = 1.225 + 0.027 X1 + 0.108 X2 + 0.014 X1X2, R 2 = 0.996, R 2 adj = 0.995 Tf Cu(DDC)2 = 1.188 − 0.127 X2, R 2 = 0.793, R 2 adj = 0.780 (8) The low magnitude of the coefficients (β1, β2, and β12) of Equations (7) and (8) suggests a minor impact of flow rate (X1) and mobile phase composition (X2) on this response. In fact, for Cu(DDC)2, the flow rate (X1) did not have a statistically significant impact. Subsequently, this response is shaped by a linear function and not characterized by a surface response model. As for DSF, the model ( Figure 6) shows a linear behavior. Overall, taking into consideration the β0 coefficient magnitude for both analytes, the criteria of this ATP specification are always met, regardless of the experimental conditions. This behavior is also reported for other compounds evaluated under different chromatographic methods [25].

Tailing Factor
Chromatographic peaks are ideally represented by Gaussian curves. However, in practice, they are usually asymmetrical due to the interaction of the analytes and the stationary phase. Thus, T f is a measure of peak symmetry [28]. According to official guidelines, and to reduce peak integration variability, T f < 2.0. For each analyte, T f can be determined by T f DSF = 1.225 + 0.027 X 1 + 0.108 X 2 + 0.014 X 1 X 2 , R 2 = 0.996, R 2 adj = 0.995 (7) T f Cu(DDC)2 = 1.188 − 0.127 X 2 , R 2 = 0.793, R 2 adj = 0.780 (8) The low magnitude of the coefficients (β 1 , β 2 , and β 12 ) of Equations (7) and (8) suggests a minor impact of flow rate (X 1 ) and mobile phase composition (X 2 ) on this response. In fact, for Cu(DDC) 2 , the flow rate (X 1 ) did not have a statistically significant impact. Subsequently, this response is shaped by a linear function and not characterized by a surface response model. As for DSF, the model ( Figure 6) shows a linear behavior. Overall, taking into consideration the β 0 coefficient magnitude for both analytes, the criteria of this ATP specification are always met, regardless of the experimental conditions. This behavior is also reported for other compounds evaluated under different chromatographic methods [25].

Critical Peak Resolution
Peak resolution assesses the degree of separation between two adjacent peaks. A resolution value of zero pinpoints a complete peak co-elution, whereas a resolution value of 1.5 is the minimum for baseline separation [28]. As per FDA recommendations, Res > 2.0. The selection of this variable arises from the visual inspection of chromatograms at 217 nm which indicates critical peak proximity between DSF and its impurities, that can be determined using: Res DSF = 3.367 − 0.15 X 1 − 4.033 X 2 + 0.15 X 1 X 2 + 0.667 X 22 2 , R 2 = 0.999, R 2 adj = 0.999 (9) 021, 9, x FOR PEER REVIEW 12 of 20

Critical Peak Resolution
Peak resolution assesses the degree of separation between two adjacent peaks. A resolution value of zero pinpoints a complete peak co-elution, whereas a resolution value of 1.5 is the minimum for baseline separation [28]. As per FDA recommendations, Res > 2.0. The selection of this variable arises from the visual inspection of chromatograms at 217 nm which indicates critical peak proximity between DSF and its impurities, that can be determined using: Res DSF = 3.367 − 0.15 X1 − 4.033 X2 + 0.15 X1X2 + 0.667 X22 2 , R 2 = 0.999, R 2 adj = 0.999 (9) According to the coefficients of Equation (9), this CAA is mainly governed by the mobile phase composition (X2). As expected, an increase in the % of ACN strongly reduces peak Res due to the lower interaction between the analytes and the stationary phase, promoted by the stronger nonpolar conditions of the mobile phase. The remaining terms are not relevant due to their reduced magnitude. Nonetheless, the model shows some degree of curvature (Figure 7). This trend is also observed for other compounds evaluated under different chromatographic conditions [25].  According to the coefficients of Equation (9), this CAA is mainly governed by the mobile phase composition (X 2 ). As expected, an increase in the % of ACN strongly reduces peak Res due to the lower interaction between the analytes and the stationary phase, promoted by the stronger nonpolar conditions of the mobile phase. The remaining terms are not relevant due to their reduced magnitude. Nonetheless, the model shows some degree of curvature (Figure 7). This trend is also observed for other compounds evaluated under different chromatographic conditions [25].

Critical Peak Resolution
Peak resolution assesses the degree of separation between two adjacent peaks. A resolution value of zero pinpoints a complete peak co-elution, whereas a resolution value of 1.5 is the minimum for baseline separation [28]. As per FDA recommendations, Res > 2.0. The selection of this variable arises from the visual inspection of chromatograms at 217 nm which indicates critical peak proximity between DSF and its impurities, that can be determined using: Res DSF = 3.367 − 0.15 X1 − 4.033 X2 + 0.15 X1X2 + 0.667 X22 2 , R 2 = 0.999, R 2 adj = 0.999 (9) According to the coefficients of Equation (9), this CAA is mainly governed by the mobile phase composition (X2). As expected, an increase in the % of ACN strongly reduces peak Res due to the lower interaction between the analytes and the stationary phase, promoted by the stronger nonpolar conditions of the mobile phase. The remaining terms are not relevant due to their reduced magnitude. Nonetheless, the model shows some degree of curvature (Figure 7). This trend is also observed for other compounds evaluated under different chromatographic conditions [25].

Retention Time
Retention time is a measure of the time each analyte spends in the stationary and in the mobile phase [28]. Generally, it should be as low as possible in order to reduce the analysis time and costs of organic solvents. For each analyte, R t can be determined by R t DSF = 4.028 − 1.036 X 1 − 2.92 X 2 + 0.61 X 1 X 2 + 1.585 X 22 2 , R 2 = 0.996, R 2 adj = 0.995 (10) R t Cu(DDC)2 = 6.95 − 2.127 X 1 − 8.557 X 2 + 1.763 X 1 X 2 + 5.115 X 22 2 , R 2 = 0.996, R 2 adj = 0.995 (11) According to Equations (10) and (11), both critical method parameters have a negative impact on this parameter, as indicated by the negative β 1 and β 2 coefficients. Although flow rate (X 1 ) has a strong influence in reducing drug R t , an increase in the organic phase (X 2 ) is the main factor that reduces the interaction of the analytes with the stationary phase, as indicated by the highest magnitude of β 2 coefficients. In addition, according to the β 22 coefficient, the surface response models display curvature (Figure 8). This behavior is also seen for other compounds and chromatographic methods [25].

Retention Time
Retention time is a measure of the time each analyte spends in the stationary and in the mobile phase [28]. Generally, it should be as low as possible in order to reduce the analysis time and costs of organic solvents. For each analyte, Rt can be determined by Rt DSF = 4.028 − 1.036 X1 − 2.92 X2 + 0.61 X1X2 + 1.585 X22 2 , R 2 = 0.996, R 2 adj = 0.995 (10) Rt Cu(DDC)2 = 6.95 − 2.127 X1 − 8.557 X2 + 1.763 X1X2 + 5.115 X22 2 , R 2 = 0.996, R 2 adj = 0.995 (11) According to Equations (10) and (11), both critical method parameters have a negative impact on this parameter, as indicated by the negative β1 and β2 coefficients. Although flow rate (X1) has a strong influence in reducing drug Rt, an increase in the organic phase (X2) is the main factor that reduces the interaction of the analytes with the stationary phase, as indicated by the highest magnitude of β2 coefficients. In addition, according to the β22 coefficient, the surface response models display curvature (Figure 8). This behavior is also seen for other compounds and chromatographic methods [25].

Capacity Factor
Although similar to the drug Rt, the capacity factor, or retention factor, represents the degree of retention of the analyte in the column and takes into consideration the void time (i.e., the retention of an unretained analyte) [28]. Usually, when using isocratic conditions, 1 < k' < 20, as higher k' values lead to extensive peak broadening. This parameter becomes particularly relevant when using increased proportions of organic solvents, as they reduce the retention of the analytes. As per FDA, k' > 2, and can be determined using: k'DSF = 2.817 − 2.762 X2 + 1.416 X22 2 , R 2 = 0.996, R 2 adj = 0.996 (12) k'Cu(DDC)2 = 5.463 − 8.067 X2 + 4.655 X22 2 , R 2 = 0.997, R 2 adj = 0.996 (13) According to Equations (12) and (13), flow rate (X1) does not have an impact on the capacity factor. In fact, k' is independent of flow rate and column length [28]. Therefore, this parameter is explained by a simpler polynomial function and not characterized by a surface response model. The negative sign of the β2 coefficient indicates that an increase in the organic phase (X2) strongly reduces the capacity factor, thus eluting the analytes in less time.

Capacity Factor
Although similar to the drug R t , the capacity factor, or retention factor, represents the degree of retention of the analyte in the column and takes into consideration the void time (i.e., the retention of an unretained analyte) [28]. Usually, when using isocratic conditions, 1 < k' < 20, as higher k' values lead to extensive peak broadening. This parameter becomes particularly relevant when using increased proportions of organic solvents, as they reduce the retention of the analytes. As per FDA, k' > 2, and can be determined using: k' DSF = 2.817 − 2.762 X 2 + 1.416 X 22 2 , R 2 = 0.996, R 2 adj = 0.996 (12) k' Cu(DDC)2 = 5.463 − 8.067 X 2 + 4.655 X 22 2 , R 2 = 0.997, R 2 adj = 0.996 (13) According to Equations (12) and (13), flow rate (X 1 ) does not have an impact on the capacity factor. In fact, k' is independent of flow rate and column length [28]. Therefore, this parameter is explained by a simpler polynomial function and not characterized by a surface response model. The negative sign of the β 2 coefficient indicates that an increase in the organic phase (X 2 ) strongly reduces the capacity factor, thus eluting the analytes in less time.

Method Development and Optimization: MODR
Aiming at obtaining a validated method, as per ICH and FDA guidelines, i.e., a chromatographic method that complies with the established requirement for the CAAs, it is necessary to determine the method operable design region (MODR). Every condition is associated with a certain level of failure risk, expressed as the probability of failure (%), calculated using Monte Carlo simulations. This risk arises not only from the method parameters and their variation but also from unexpected and random sources. As such, a random standard deviation of 0.05 and 0.5 for X 1 and X 2 , respectively, was assumed. The identification of this region is depicted in white, in Figure 9, based on the criteria defined in Table 5. Each point of the white region of the MODR represents a possible chromatographic Chemosensors 2021, 9,172 14 of 20 method (in terms of flow rate (X 1 ) and mobile phase composition (X 2 )), compliant with the ATP (Table 1), and with a risk of failure below 5%. necessary to determine the method operable design region (MODR). Every condition is associated with a certain level of failure risk, expressed as the probability of failure (%), calculated using Monte Carlo simulations. This risk arises not only from the method parameters and their variation but also from unexpected and random sources. As such, a random standard deviation of 0.05 and 0.5 for X1 and X2, respectively, was assumed. The identification of this region is depicted in white, in Figure 9, based on the criteria defined in Table 5. Each point of the white region of the MODR represents a possible chromatographic method (in terms of flow rate (X1) and mobile phase composition (X2)), compliant with the ATP (Table 1), and with a risk of failure below 5%. In order to validate the predictability of the operable region, four random conditions ( Figure 9, green dots) were experimentally evaluated and compared to the predicted CAA responses by the models (Table 5). For these levels, the MODR also complies with the ATP (Table 1), supporting the good predictability of the models. For convenience, the chosen method parameters to proceed for formal validation were set as being 70:30 (v/v) ACN:H2O, at a flow rate of 1.0 mL/min (Figure 9, blue dot). These conditions guarantee the suitability of the chromatographic method and are potentially associated with a better ruggedness, i.e., present a lower probability of failure when used at another laboratory. Moreover, regions close to the limits of X1 or X2 should be avoided, as they are on the edge of the knowledge space.  In order to validate the predictability of the operable region, four random conditions (Figure 9, green dots) were experimentally evaluated and compared to the predicted CAA responses by the models (Table 5). For these levels, the MODR also complies with the ATP (Table 1), supporting the good predictability of the models. For convenience, the chosen method parameters to proceed for formal validation were set as being 70:30 (v/v) ACN:H 2 O, at a flow rate of 1.0 mL/min (Figure 9, blue dot). These conditions guarantee the suitability of the chromatographic method and are potentially associated with a better ruggedness, i.e., present a lower probability of failure when used at another laboratory. Moreover, regions close to the limits of X 1 or X 2 should be avoided, as they are on the edge of the knowledge space.

System Suitability
System suitability testing aims at assessing if the chromatogram and the HPLC equipment are able to produce accurate and precise results [21]. For that, several parameters, including the previously identified CAAs, are evaluated. The precision of the retention times and peak areas were compliant for both analytes, according to the RSD ≤ 2.0% for the three quality control standards, thus supporting the system's ability for their detection and quantification (Table 6). Furthermore, the efficiency of the column is also evidenced for all parameters. According to the ATP of the established AQbD strategy (Table 1), compliance with these CAAs was a pre-validation requisite.

Detection and Quantification Limits
The lowest concentration at which each analyte can be detected or identified with good accuracy and precision was determined according to the standard deviation of the response and slope of the linear regression of the calibration curve. The estimated LOD for DSF and Cu(DDC) 2 were 0.034 and 0.028 µg/mL, respectively. Additionally, the LOQ for DSF and Cu(DDC) 2 were 0.10 and 0.09 µg/mL, respectively.

Linearity
Linearity of the calibration standards was observed over the range of 0.1 to 100 µg/mL, according to the regression equations and correlation coefficients obtained through the least squares method.
The mean equations for both analytes are presented below (Equations (14) and (15), n = 6). A 1/x 2 weighting factor was employed due to the heteroscedasticity of the samples.

Accuracy and Precision
The evaluation of the accuracy and precision of the method is a further guarantee of the quality of the calibration curves and should be observed at intra-and interday levels. According to Table 7, both parameters are verified since the deviations to nominal concentration and among independent samples did not exceed 15. As such, these results corroborate the accuracy, reliability, and reproducibility of the method.

Robustness and Ruggedness
The assessment of the robustness and ruggedness of chromatographic methods aims at ensuring the reliability of the results. To study the method's robustness, the ICH suggests the evaluation of small method variations, such as modifications to the mobile phase composition, pH, different columns, oven temperature, and flow rate [19,29]. According to the risk assessment ( Figure 4 and Table 3), mobile phase composition and flow rate were identified as critical process parameters due to their potential high criticality. The systematic study, using Design of Experiments (Section 2.5.3) and the wide MODR region, support the robustness of the method. In order to study the method's ruggedness, the FDA takes into consideration the expected variations occurring between different laboratories. Consequently, the optimized method was evaluated in a different apparatus and column. Despite the expected minor variations, the different conditions met the CAA criteria for system suitability, including a good critical peak resolution, capacity factor, theoretical plate number, and tailing factor (Table 8). Additionally, the method complies with the requirements for repeatability (intraday precision) and accuracy (Table 9). 3.4.6. Specificity The specificity of the method is a parameter that ensures the correct determination of the analytes in the presence of possible interferents. For that, the supernatant of blank nanoparticles was analyzed using the optimized chromatographic method conditions and compared to a representative chromatogram of a 0.2 µg/mL quality control solution (Figure 10a). No peaks were observed at the retention times of the analytes, indicating the method can be regarded as specific. 3.4.6. Specificity The specificity of the method is a parameter that ensures the correct determination of the analytes in the presence of possible interferents. For that, the supernatant of blank nanoparticles was analyzed using the optimized chromatographic method conditions and compared to a representative chromatogram of a 0.2 µg/mL quality control solution (Figure 10a). No peaks were observed at the retention times of the analytes, indicating the method can be regarded as specific.

Stability
Stability testing, as part of the chromatographic method validation, aims at ensuring there is no degradation or adhesion to glassware during sampling procedures and laboratory analysis. As such, stability studies were performed at two levels: stock solution

Stability
Stability testing, as part of the chromatographic method validation, aims at ensuring there is no degradation or adhesion to glassware during sampling procedures and laboratory analysis. As such, stability studies were performed at two levels: stock solution (freeze-thaw and intermediate-term stability) and quality control solutions testing (autosampler and 72 h stability, at 4 • C). According to the results (Table 10), both analytes met the established accuracy and precision criteria in these conditions. Consequently, the quality of the data is guaranteed.

Method Applicability: Nanostructured Lipid Carriers
The optimized method was applied to determine the content of DSF and Cu(DDC) 2 in lipid nanoparticles, with a mean particle size of 140 ± 4 nm, a polydispersity index of 0.131, and a zeta potential of −40 mV ± 2 mV. Overall, the entrapment efficiency values of DSF and Cu(DDC) 2 were found to be 99.7 ± 0.2% and ca. 100%, as no peak was identified in the aqueous phase of the NLCs corresponding to Cu(DDC) 2 (Figure 10b). The determined loadings of DSF and Cu(DDC) 2 were 2.10 ± 0.04% and 6.6 ± 0.2 ‰. Despite the reduced DL of Cu(DDC) 2 , this compound has been described as highly potent at decreasing cancer cell viability, displaying activity in nanomolar concentrations against glioblastoma, breast and lung cancer [4]. This method may also be used for other routine tests that include the use of these analytes.

Conclusions
An AQbD strategy was implemented for the simultaneous quantification of disulfiram and its anticancer compound, copper (II) diethyldithiocarbamate, by liquid chromatography. According to risk assessment, the method parameters with a high impact on method performance (by means of the theoretical plate number, critical resolution, retention time, tailing factor, and capacity factor) were evaluated. Under the umbrella of AQbD, the application of a face-centered central composite design allowed the characterization of these parameters. As such, the multidimensional region that includes the method parameters that best describe the suitability of the chromatographic method (with a negligible risk of failure) was defined using Monte Carlo simulations and experimentally validated. The ideal conditions (a mobile phase containing acetonitrile:water (70:30, v/v) at 1 mL/min) were validated according to regulatory guidelines in terms of linearity, accuracy, precision, specificity, robustness, and stability. In addition, the method was evaluated for ruggedness, using different HPLC equipment and column, meeting the criteria for system suitability, intraday accuracy, and precision (repeatability). Overall, the validation of the method supports the implementation of AQbD in defining chromatographic conditions for research or industrial applications, exemplified here with the quantification of the analytes in nanostructured lipid carriers.