Improving image contrast and material discrimination with nonlinear response in bimodal atomic force microscopy

Atomic force microscopy has recently been extented to bimodal operation, where increased image contrast is achieved through excitation and measurement of two cantilever eigenmodes. This enhanced material contrast is advantageous in analysis of complex heterogeneous materials with phase separation on the micro or nanometre scale. Here we show that much greater image contrast results from analysis of nonlinear response to the bimodal drive, at harmonics and mixing frequencies. The amplitude and phase of up to 17 frequencies are simultaneously measured in a single scan. Using a machine-learning algorithm we demonstrate almost threefold improvement in the ability to separate material components of a polymer blend when including this nonlinear response. Beyond the statistical analysis performed here, analysis of nonlinear response could be used to obtain quantitative material properties at high speeds and with enhanced resolution.

D etailed analysis of cantilever dynamics in atomic force microscopy (AFM) allows one to go beyond simple topography mapping to study the physical properties of a material's surface 1,2 . The recent development of multifrequency AFM facilitates such analysis by measuring response at harmonics 3,4 , in a continuous frequency band 5 , at discrete tones near one resonance 6 and at two or more flexural resonances [7][8][9] . This later so-called bimodal AFM has demonstrated increased material contrast 10 and material property mapping 11 , but the limited number of measured signals restricts mapping to parameters of very simple tipsurface interaction models. Additional signals are available if one measures response at the intermodulation products or mixing frequencies of the bimodal drive [12][13][14] but the utility of these extra signals for imaging has never been demonstrated. We show that image contrast at mixing frequencies can be larger than at drive frequencies, and we use Fisher's linear discriminant analysis 15 on the collective multifrequency response to demonstrate quantitative improvement of material contrast.
Traditional modes of dynamic AFM derive each image point from two measured quantities: the amplitude and phase of cantilever oscillation at the frequency of the drive. Feedback adjusts the probe height so as to keep the oscillation amplitude constant and this constant amplitude height image is typically interpreted as surface topography. Changes in the phase indicate a variation of material property. Assuming sinusoidal motion, the phase image can be interpreted as a map of energy dissipation 16 and for many years this phase image was the only information the AFM operator had to understand the physical character of the surface.
To improve the AFM's ability to sense material properties one must increase the number of measured signals. Multifrequency AFM 17 approaches this task in the frequency domain by exciting and measuring the response of the cantilever at many frequencies in the time required to record one pixel, where each frequency provides two observable quantities. Bimodal AFM 7 , where the cantilever is excited at the resonant frequencies of two different eigenmodes, has been extended to multiple eigenmodes 9,18,19 . This approach gives more signals, but it requires broadband detection. With the photodetectors used in today's AFMs it is typically not feasible to measure more than the first few modes.
Rather than extending the number of eigenmodes, more information can be obtained in the limited detection bandwidth by capturing response at non-driven frequencies. When the cantilever is driven with two tones at frequencies f 1 and f 2 , the nonlinear tip-surface interaction generates harmonics and so-called intermodulation or mixing products at frequencies where n and m are integers (positive or negative), and |n| þ |m| is the order of the intermodulation product 6 . Below we show that this nonlinear response in bimodal AFM can be utilized to produce higher image contrast than the linear response.

Results
Free and engaged cantilever dynamics. An AFM cantilever was inertially excited using a shaker piezo driven with a superposition of two pure tones at f 1 ¼ 78.5 and f 2 ¼ 500.5 kHz, very close to the first two flexural eigenfrequencies. The excitation was not exactly on resonance, as the second resonance frequency is not a rational fraction of the first. True bimodal oscillation is therefore incommensurate, with each oscillation at the first flexural eigenfrequency slightly different from the previous. To ensure periodicity in the measurement window or pixel time T, the two drive frequencies were chosen to be integer multiples of one base frequency, Df ¼ 1/T ¼ 0.5 kHz. This choice avoids spectral leakage and enables the measurement of both amplitude and phase at intermodulation frequencies. When free from the surface, the drive signal was adjusted such that the amplitude of tip motion at the first eigenmode was 17 nm (290 mV). The free amplitude at the second eigenmode was adjusted to be 20% of the voltage amplitude at the first mode, as measured by the photodetector. Using the theoretical difference in detector responsivity, we estimate the free oscillation amplitude of tip motion at the second mode to be about 1 nm (see Methods for calibration details). The spectrum of the freely oscillating cantilever showed weak mixing products and harmonics (see Fig. 1a), which we attribute to nonlinearity in either the drive electronics, the excitation piezo, the photodetector or the digital sampling electronics. This 'background' nonlinearity can limit the ability to precisely reconstruct tip-surface force, but it will not effect the statistical analysis used here, which does not require knowledge of the physical theory connecting tip-surface force with cantilever motion. When the oscillating cantilever engaged a surface, many new mixing products appeared above the background (compare Fig. 1a,b) clearly due to the tip-surface interaction. The response at 6f 1 is stronger than other harmonics as this frequency  Table 1 plotted together with the approximate response function of each eigenmode. Although the force responsivity and PD amplitude are much higher on resonance, off resonance contrast is excellent, in many cases better than on resonance.
falls closest to the second eigenfrequency of the cantilever. The difference frequency is a factor in many of the stronger intermodulation products (see Table 1).
Imaging. We measured the response at the dominant peaks in the spectra, listed in Table 1, while scanning the surface of a blend of polystyrene and low-density polyethylene (PS-LDPE), Fig. 2.
A multifrequency lockin amplifier 20 captured the amplitude and phase at all frequencies in real time during a single scan. As LDPE is much softer than PS (B0.1 GPa versus B2 GPa 21 ), excellent contrast is seen in the amplitude and phase images at all frequencies, except the amplitude at f 1 used for feedback (see Supplementary Table I).
To quantify image contrast at different frequencies we converted the amplitude and phase into their real and imaginary quadratures ( Supplementary Fig. 1) and made image histograms of these quantities ( Supplementary Fig. 2). Representation of the response data in a Cartesian plane avoids ambiguity associated with the wrapping of phase, and it allowed all histograms to have the same physical units of detector voltage. The image histograms contained one peak corresponding to PS and one peak for LDPE. We fit the histograms to a binormal distribution (see Methods) and define a contrast metric, Where, m and s 2 are the mean and variance of the two normal distributions. Two contrast values c real and c imag were calculated at each frequency, from the real and imaginary quadrature images, respectively. Although the response amplitudes at the driven frequencies f 1 and f 2 were orders of magnitude larger than at harmonics and mixing tones, the latter often showed higher contrast (see Fig. 1c and Table 1). In other words, discrimination of the two polymer components was greater at weakly responding mixing frequencies than at strongly responding drive frequencies.
A second sample consisting of a blend of polystyrene and polymethylmethacrylate (PS-PMMA) showed lower contrast than PS-LDPE at all frequencies (see Supplementary Figs 3,4  and Supplementary Tables II,III), which was expected as the two polymers have similar elastic moduli 22 . For this sample the largest contrast was found at f 1 , with several mixing frequencies having higher contrast than that at f 2 .
Linear discriminant analysis. One scan generates a 34-dimensional data set consisting of the complex amplitudes at 17 frequencies.
One can interpret such high-dimensional data sets with physical models or black-box models. Physical modelling uses analytic or numeric methods to understand the data in terms of underlying physical laws. Normal bimodal AFM, with its 4-dimensional data set (two frequencies) restricts physical interpretation to simple models 11 . The lack of a good method to calibrate the cantilever and measure the actual tip motion at the second eigenmode also limits the use of physical modelling in bimodal AFM 23 . In contrast, black-box modelling does not rely on a physical theory of the system being imaged, but rather attempts to draw conclusions from correlations in the raw data. Principle component analysis and neural networks have been applied to multidimensional AFM data 24,25 . We use Fisher's linear discriminant analysis (LDA) 15 on two polymer components to demonstrate a quantitative improvement in image contrast by including response at harmonics and intermodulation products. Each pixel is an independent measurement of 34 'features' represented by a vector in the 34-dimensional space. We require that each pixel belong to one of two classes, for example, PS or LDPE. The LDA algorithm will find the linear combination of the measurement bases, which maximizes the separation between the two classes and minimizes the variance within the classes 26 . This optimal projection is calculated from manually selected training data having an assumed classification. A new image is formed by applying the optimal projection to all pixels (Fig. 3a,b,d,e). Three linear discriminant analyses were performed using image features defined in three different ways: one using response only at the two driven frequencies (denoted 'dual'); one using all 17 frequencies ('all'); and one using only the first drive frequency f 1 ('single'). This last analysis allows for comparison of bimodal AFM with standard single frequency AM-AFM. A separate scan performed with the second drive frequency turned off gave the same result as analysing only f 1 with a bimodal drive.
Histograms for the projected images were fit to a binormal distribution and the contrast metric equation (3) was calculated (Fig. 3c,f). A small increase of contrast was found going from single frequency to bimodal AFM, 13% increase for the PS-LDPE and 5% increase for PS-PMMA. However, when including all frequencies in the analysis a much larger improvement in contrast was found: almost threefold for PS-LDPE and 15% increase for PS-PMMA. Similar improvements were found in multiple scans with several different cantilevers, including a much stiffer and higher Q-factor cantilever (BudgetSensor Tap300-Al, k ¼ 32 N m À 1 , Q 1 ¼ 548). These improvements were obtained from measuring response at additional non-driven frequencies, without changing experimental conditions such as the driving force or feedback set-point. Therefore, contrast was gained without increasing the applied force to the surface.

Discussion
The LDA analysis shows that including high order intermodulation products enhances the ability to discriminate between the two components of the polymer blend. The image contrast is influenced by at least two factors: the physical process creating the response and noise when measuring the response. The cantilever is more sensitive to force on resonance, resulting in a higher signal-to-noise (SNR) ratio. However, the physical process that gives rise to contrast is the tip-surface force, which is a nonlinear perturbation of the cantilever dynamics. This nonlinear process can be much stronger off resonance, at harmonic and intermodulation frequencies.
Combining nonlinearity and resonance is most favourable. For instance, more mixing terms with larger SNR would occur if the second eigenmode resonance frequency was adjusted such that D ¼ f 2 À 6f 1 was smaller. However, if the second resonance happened to be at exactly six times the first, it is not clear that choosing f 2 ¼ 6f 1 would be advantageous. Driving at the 6th harmonic of f 1 would not create mixing frequencies, only harmonics (integer multiples of f 1 ). In this case it would be advantageous to drive the cantilever at f 1 and f 2 ¼ 6f 1 ± Df, where Df correspond to the measurement bandwidth for one pixel, such that many intermodulation frequencies occur close to resonance 6 .
Thus new cantilever designs could be used to enhance cantilever response 27 and low noise detectors would help improve the SNR off resonance [28][29][30] . The enhanced contrast clearly seen at several mixing frequencies points towards further possible development of the methods used here. The use of blackbox models such as LDA can help digest the high-dimensional data sets obtained with this and other emerging multifrequency AFM methods. For complex samples LDA extends to more than two classes and if the use of training data is not possible, there exists a large variety of unsupervised clustering algorithms, which could be applicable 31 . If calibration methods for higher modes are developed, physical models can be formulated, which take into account frequency mixing with multiple eigenmodes. Better quantitative analysis of the forces on the sample and therefore better discrimination of material composition and changes in topography will be possible. We recently described a method starting from an arbitrary physical model to approximate material properties (model parameters) using response at mixing frequencies 32 and multiple eigenmodes 33 .

Methods
Sample preparation. The PS-LDPE sample was obtain from Bruker Corporation (HarmoniX test sample) and the PS-PMMA sample was made from polymers obtained from Sigma Aldrich. 17 mg PS (M w ¼ 280,000 kDa) and 52 mg PMMA (M w ¼ 120,000 kDa) were dissolved in 8 ml toluene and stirred for 4 h. A silicon substrate with native oxide was sonicated first in acetone, then isopropanol, followed by 4 min of O 2 plasma. The polymer solution was spin cast on the substrate at 500 rpm for 5 s and 3000 rpm for 60 s, forming a B20-nm thick layer.   AFM imaging. AFM measurments were performed on a Nanowizard 3 BioScience AFM (JPK Instruments, Germany) with an additional multifrequency lockin amplifier 20 (Intermodulation Products AB, Sweden). Thermal noise power spectra were used to determine the two first flexural eigenmode resonance frequencies and quality factors of a RFESPA AFM cantilever (Bruker):f 1 ¼ 78:973 kHz, f 2 ¼ 501:055 kHz and Q 1 ¼ 183, Q 2 ¼ 485 respectively. The first mode was calibrated using the Sader-Higgins method 34,35 to obtain a dynamic cantilever stiffness k 1 ¼ 2.5 N m À 1 and the inverse photodetector responsivity a 1 À 1 ¼ 58 nm V À 1 . For the second mode ideal results from Euler-Bernoulli beam theory were used 17 to estimate k 2 ¼ 39.32k 1 ¼ 99 N m À 1 and a 2 À 1 ¼ a 1 À 1 / 3.573 ¼ 17 nm V À 1 . Using this method the estimated quality factor of the second mode was a factor of 3 larger than the measured Q 2 . This discrepancy highlights difficulties in calibration of the second mode.
To avoid Fourier leakage, the drive frequencies f 1 ¼ 78.5 and f 2 ¼ 500.5 kHz were chosen as integer multiples of a common base frequency Df ¼ 500 Hz. Df constituted the measurement bandwidth or inverse of the measurement time window T ¼ 1/Df. The discrete Fourier transform or a lockin calculation (Fourier sum for one frequency) was performed over this time window. For T ¼ 2 ms per pixel and 256 pixels per line, the scan rate is B1 line per second, typical for dynamic AFM imaging. Feedback was used to adjust the surface position so as to keep the amplitude at f 1 constant Numerical methods. All histograms were normalized (total area equalled unity) and a sum of two normal distributions was fit to the histograms, where m 1,2 and s 2 1;2 are the mean and variances for the two distributions 1 and 2, p is the ratio of pixels in the first distribution and (1 À p) the ratio of pixels in the second. The contrast metric was calculated from (3). Fisher's linear discriminant analysis was calculated using ref. 26 equation (4.30).