One factor at a time is a procedure with an assumption buried inside it. You hold the method fixed, walk one variable across a few settings, keep the value that gave the best chromatogram, lock it, and move to the next variable. That sequence recovers the true optimum only if the response is separable — if the best column temperature is the best column temperature whatever the ion-pairing agent is doing, and the best ion-pairing concentration is the same at every pH. In a peptide reversed-phase gradient it is not.

The ion-pairing factor demonstrates this on its own. Working with 18-residue amphipathic α-helical model peptides, raising trifluoroacetic acid concentration or adding sodium perchlorate to a phosphoric acid mobile phase at pH 2.0 improved peak shape and increased retention through ion-pairing with the positively charged peptides; at pH 7.0 the addition of salt did little beyond improving peak shape [1]. The ion-pairing factor is therefore large or negligible depending on which pH you happened to fix first. An OFAT walk that settled pH at 7.0 early would report, correctly for the conditions it tested, that adding salt does little beyond improving peak shape.

Interactions do more than mis-size an effect; they reorder peaks. In a supercritical fluid chromatography separation of terbutaline enantiomers and their monosulfate conjugates, removing the mobile phase additive switched the elution order of the two peak pairs, and the response surface favoured high temperature and high flow rate — up to 2.25 mL/min — compared with what previous reports had indicated. The authors state that these multifactorial interactions would not have been detected by the common one-factor-at-a-time approach [2].

The second failure is geometric. OFAT climbs along the coordinate axes, so where the ridge of good performance runs diagonally across two factors, an axis walk stops short and reports a plateau. Peptide ion-pairing supplies a worked example of an interior optimum that convention had missed. TFA is generally used at 0.05–0.1% (6.5–13 mM). Re-examined on three groups of 10-residue synthetic peptides carrying one, three or five positive charges, at 25 and 70 °C, that traditional range was found not to be optimum for many peptide applications: 0.2–0.25% TFA gave optimum resolution for mixtures containing multiply charged peptides, and high temperature was an effective complement to those TFA levels [3]. Both halves matter: the optimum sits inside the factor range rather than at its edge, and the same study names a second factor, high temperature, as an effective complement to those TFA levels.

Two stages asking two different questions

Screening and optimisation are not the same experiment run at different resolutions. Screening asks which of a long candidate list moves the response at all, and is sized to answer that and nothing more. Optimisation asks what the surface looks like over the survivors, and needs enough distinct levels per factor to fit curvature.

The design families sort accordingly. Plackett–Burman designs are widely used in validation studies to show the absence of significant effects; fractional factorial designs and their extensions such as central composite designs are the most popular optimisers; Box–Behnken and Doehlert designs are increasingly used as efficient alternatives; and where the factor combinations a design calls for cannot be realised in practice, or the total number of runs is capped, D-optimal designs can be very powerful [4].

Run counts make the division concrete. In one RP-HPLC optimisation, candidate method attributes were screened by a Pareto-chart-based fractional factorial design, and the survivors were carried into a Box–Behnken design of 17 runs; analysis of variance, lack-of-fit and response-surface plots then evaluated retention time, tailing factor and theoretical plate count [5]. In a UHPLC separation of sixteen compounds, screening designs identified the critical parameters first — in ethanol-based mobile phases temperature turned out to be critical for both retention and separation — and response surface methodology was run only on those [6].

The factors in a peptide gradient, and why they are coupled

For peptides the shortlist is short and well characterised: gradient slope, column temperature, ion-pairing agent concentration, mobile phase pH and the stationary phase. A study spanning 165 peptides totalling 1698 amino acid residues varied column temperature, gradient slope, ion-pairing reagent concentration and stationary phase type, and converted the retention times into per-residue retention coefficients. TFA concentration had the strongest impact on the retention of peptides containing increasing numbers of basic residues — arginine, lysine, histidine — while coefficients for other residues were minimally affected. Raising the separation temperature moderately decreased retention coefficients, with the exception of isoleucine, leucine, valine and proline, which the authors suggest enhance peptide retention at elevated temperature [7].

That is the mechanism behind the coupling. The ion-pairing factor acts through peptide charge; the temperature factor acts partly through conformation, and therefore acts differently on different sequences. Across 5–80 °C, amphipathic α-helical peptide analogues showed a greater effect of temperature on elution behaviour than random-coil analogues, attributed to unfolding of the α-helix as temperature rose during the separation [8]. Two peptides in the same injection can have retention-versus-temperature relationships of different shape, which is why a temperature coefficient measured on one peak does not carry to the next.

Pairing gradient slope with ion-pairing concentration has already been used as a two-factor design in peptide work. Adsorption isotherms for three peptides were determined at eight TFA concentrations from 2.6 to 37.3 mM (0.02–0.29 vol-%) on an end-capped C18 phase, with fluorine–proton nuclear Overhauser experiments providing direct evidence of ion-pair formation between a peptide and TFA in acetonitrile–water. At low TFA the elution profiles were distorted, showing a spike at the front and a long tailing rear, because the free peptide ion has strong secondary interactions and very low saturation capacity; the practical method-development discussion rests on an experimental design using gradient slope and TFA concentration as factors [9]. Scope: those isotherms were measured under overloaded conditions, so the peak-shape result speaks most directly to loading-limited work.

What the polynomial is fitting

A response-surface design fits an empirical polynomial, usually second order, and the second-order term earns its place because the first-order one demonstrably fails. In reversed-phase separation of monoclonal antibodies and their fragments, gradient steepness and temperature could not be optimised using van 't Hoff type linear models, and the common linear solvent strength model also generated some error in predicting retention times; quadratic models brought the relative error on predicted retention time to between 0.5 and 1%, from nine experiments for the gradient-time-and-temperature optimisation, or six, taking six to eight hours [10]. Those figures are for intact mAbs and fragments on wide-pore phases at 60–90 °C, not for small synthetic peptides, but the geometry of the problem is the same one a peptide method meets.

Curvature in the temperature axis is not occasional. Across 1211 van 't Hoff plots, each a solute–mobile phase–day combination drawn from more than 25 000 retention factors measured on 101 small molecules, 98% were fitted significantly better by an expression carrying a quadratic term in inverse temperature than by the linear form [11]. A two-level design in temperature cannot see any of that.

There is also a response that polynomials handle badly, and it is the one most people reach for first. Resolution is not a smooth function of the factors when peaks cross. Applied over 1–60 mM of phosphoric acid, TFA, pentafluoropropionic acid and heptafluorobutyric acid to two mixtures of 18-residue synthetic peptides, one of uniform +4 net charge and one of varying +1 to +4 net charge, all four reagents increased peptide retention time with increasing concentration, to extents ordered by the hydrophobicity of the anion; above 10 mM HFBA the elution order of the four peptides reversed to +1 < +2 < +3 < +4, the opposite of the order the other three reagents produced across their entire concentration range [12]. A quadratic in resolution fitted across a region containing that crossover is fitting a fold rather than a surface: resolution passes through zero at coalescence and rises again beyond it, and a smooth quadratic fitted across it will round the crossing off rather than reproduce it.

The repair is to model something monotonic and derive resolution afterwards. Faced with elution-order changes across the experimental domain, one UHPLC development modelled the time interval between every possible pair of peaks over the whole three-dimensional domain and built the desirability function on the smallest predicted interval, arriving at a gradient slope of 2.57%/min, pH 4.85 and 33.7 °C, with resolutions above 2.4 and robustness validated experimentally [6].

Replication, randomisation and a lack-of-fit test that means something

A lack-of-fit test compares the variance the model leaves unexplained against the variance the measurement produces on its own. Without genuine replicates there is no estimate of the second term, so the ratio has no denominator and the test is simply unavailable. A replicate estimates whatever it was allowed to vary: repeat injections from one prepared vial and one column equilibration bound a narrower quantity than repeat runs set up independently at the same nominal factor settings.

How much pure error there is, and how it grows with elapsed time, is measurable. In the 101-compound retention-model evaluation above, pure experimental uncertainty in the retention factor was approximately 0.22% for measurements made on the same day, and 2.1% for measurements made on different days [11]. Close to an order of magnitude separates them, and that gap is the argument for randomising run order. Any factor whose levels happen to be walked slowly, in blocks that align with days, carries the day-to-day component in its estimated effect, and replication within a day will not expose it. Randomisation does not remove drift; it converts drift from a bias loaded onto one factor into noise spread across all of them.

The same study makes a sharper point about what a good fit certifies. Three published retention equations were fitted to the data by least squares and compared using lack-of-fit residuals, and the adsorption-based model gave the best overall fits. Yet there was virtually no correlation between a model's overall lack of fit to the experimental ln(k) and the corresponding lack of fit of the linear van 't Hoff expression, from which the authors conclude that temperature-dependent enthalpy is apparently not the cause of a model's ability to fit ln(k) [11]. Read that as a warning: a model can fit the response you measured without its extra flexibility corresponding to the mechanism you assumed it was capturing.

Design space: what is modelled, what is proven

The concept is a multidimensional region of factor settings within which the method meets its criteria. ICH Q8 presents the design space of a process as the set of factor settings providing satisfactory results, but proposes no practical methodology to define, derive or compute it [13]. Two computational routes fill that gap; both return a probability, not a guarantee.

The Bayesian route combines design of experiments with Bayesian standard multivariate regression to obtain a predictive distribution for a new response vector, and extends it to the multicriteria setting so that the predictive probability of several critical quality attributes being jointly achieved in future use of the method can be estimated [13]. The error-propagation route reaches a similar object by simulation. In an HPLC screening method developed for a 19-component mixture, gradient time, temperature and pH were the factors of a full factorial design; retention times at the beginning, apex and end of each peak were modelled by stepwise multiple linear regression, and the design space was computed as the multidimensional subspace where the probability that the separation criterion lay within acceptance limits exceeded a defined quality level, the prediction error being propagated from the modelled responses to the quality criterion by Monte Carlo simulation [14].

A modelled design space says: given this model, fitted to these runs, the probability of passing at this setting is at least p. It is conditional on the model. A proven acceptable range is a different kind of claim — you set the condition, ran it and measured the outcome — and it is correspondingly narrower, because it covers only settings actually visited.

Verification is the bridge between them, and it verifies points, not regions. A retention-modelling study across nine stationary phases performed 12 input calibration experiments per column, systematically changing gradient time, temperature and the ternary composition of the mobile phase, then overlaid the nine models to identify the parameter combinations giving at least baseline separation on all of them — a common method operational design region within which the columns were interchangeable in resolving capability. A single column-independent working point inside that region was run for verification: the experimental chromatograms agreed with the model and every column returned critical resolution above 2.0, although the elution window differed between some columns [15]. That run converts one setpoint from predicted to demonstrated and leaves the rest of the region predicted.

The same architecture has been run end to end on a peptide analyte. A UHPLC-PDA method for octreotide used risk assessment to set the critical method parameters entering a screening design, response surface methodology to model the critical method attributes over the survivors, a desirability function to locate the optimum, and Monte Carlo simulation with capability analysis to establish a robust zone within the method operable design region [16].

The model that fits its own points and misses between them

This is the characteristic failure, and it is easy to miss: every statistic on the report looks healthy. A model carrying nearly as many coefficients as there are distinct design points reproduces those points because it has enough freedom to do so, and the residuals are small because little is left over to be a residual. R² measures agreement at the runs you performed. A lack-of-fit test, where replicates make it available, tests the model against pure error at those same settings. Neither is evidence about a condition you did not run.

Three consequences follow. Add confirmation runs at interior conditions that are not design points, particularly midway between levels where a quadratic has least support, and compare predicted against measured; the verification run described above is exactly that check performed once. Do not extrapolate past the studied ranges, because the coefficients were estimated inside a box and carry no information about its outside. And treat a discontinuity in the response as a modelling error rather than a data problem: where a factor change reverses elution order, as HFBA concentration does above 10 mM [12], no smooth polynomial in resolution over that region is correct, and the fix is to change the modelled response rather than to add terms.

What the design licenses

A screening design licenses one conclusion: this factor is or is not worth carrying into the optimisation stage, over the ranges screened. It does not deliver the magnitude of that factor's effect at other settings of the other factors, because it was not built to estimate them.

A response-surface optimum licenses a setpoint to verify. Until a run at that setpoint returns the predicted chromatogram, the optimum is an inference from a fitted surface.

A modelled design space licenses a predicted probability of meeting criteria under the model, across the region the model covers. A proven acceptable range licenses the settings you ran, and nothing between them. Both are worth having, and keeping them distinct is what stops a method being run routinely at a condition that sat inside the modelled region and was never measured.

The ranges that survive this work are the ranges the system suitability limits and the method transfer protocol should be built from. A robustness range established by design, with its pure error known, states how far a receiving laboratory's temperature or ion-pairing concentration can drift before the separation becomes a different separation — which is the question a transfer is asking.