SUPPLEMENTARY INFORMATION

Sub-nanometer resolved chemical imaging by multivariate analysis of tip-enhanced Raman maps

Song Jiang, Xianbiao Zhang, Yao Zhang, Chunrui Hu, Rui Zhang, Yang Zhang, Yuan Liao, Zachary J. Smith*, Zhenchao Dong*, J.G. Hou

Content of Sections

S1. Supplementary methods: experimental setup

S2. Single-peak analysis of TERS images on molecular domains

S3. Data pre-processing methods for multivariate analysis

S4. Brief description of vertex component analysis (VCA)

S5. Estimation of the spatial resolution of TERS imaging

S6. Simulation of TERS-VCA images for a ZnTPP domain on Ag(111)

S7. Single-peak analysis of TERS images on a molecular chain

S8. Hierarchical clustering analysis (HCA) of TERS images on a molecular chain

S9. References

S1. Supplementary methods: experimental setup

Supplementary Figure 1 | Schematic drawing of our home-built experimental setup. The setup is composed of four sub-systems: a laser source for Raman excitation, a darkbox for optical filtering and alignment, a low-temperature (LT) ultrahigh-vacuum (UHV) scanning tunneling microscope (STM) for sample preparation and characterization with a built-in lens for both Raman excitation and collection, and a spectrometer equipped with a highly sensitive CCD detector for Raman spectral measurements. BS: non-polarizing beam splitter cube. ND filter: round continuously variable metallic neutral density filter.

Our TERS experiments were performed with a custom low-temperature ultrahigh-vacuum (UHV) scanning tunneling microscope (STM) under a base pressure of ~1 × 10−10 Torr at ~80 K, equipped with a side-illumination confocal optical system, as shown in FigureS1. A single-longitudinal-mode diode-pumped laser at 532 nm (CrystaLaser, CL532-100-SO) was used to provide a linearly polarized laser beam (>300:1) for Raman excitation. The beam was fiber-coupled to a dark-box via a single-mode polarization-maintaining fiber with a collimated output beam of ~1.0 mm, and then passed through a half waveplate to achieve desired p-polarization. Two reflective mirrors were used to provide freedom for optical alignment. The laser beam was introduced into the UHV chamber via a 50:50 beam splitter positioned outside the quartz viewport, which also help monitor the optical intensity of the incident beam with a power meter. A flippable 50:50 non-polarizing beam splitter was used tomonitor the focusing of the laser beam into the tunnel junction with a video camera.

The collimated beam was re-focused by an aspheric lens (fb=19.5 mm, NA=0.51) into the tip−substrate junction with an angle of 60° from the surface normal. The diameter of the laser spot on the sample surface was about 30 μm. The Raman scattered light was collected by the same lens, filtered by aRayleigh edge filter (Semrock, ultrasteep long-wave pass) to remove residual laser light, and finally fiber-coupled to a spectrometer using a slit size of 100 μm. The Raman signal light was dispersed by a 600 grooves/mm gratingand detected by a liquid-nitrogen cooled charge-coupled-device (CCD). The laser power used in our TERS experiments was about0.7 mW corresponding to a photon flux of ~100 W/cm2. (Such a small photon flux into the junction area help toavoid otherwise serious problems of molecular diffusion, desorption or even damage under laser illumination.) The spectral resolution was about 18 cm−1, limited by the instrumental resolution to maximize the collection and detection of Raman photons with a reasonable signal-to-noise ratio. All Raman spectra presented here were not corrected for the wavelength dependent sensitivity of photon detection systems.

S2. Single-peak analysis of TERS images on molecular domains

Supplementary Figure 2 | Single-peak analysis of TERS images on molecular domains.a, STM images before (left panel) and after (right panel) TERS mapping on molecular domains (–1 V, 5 pA). b,Line profile analysis for the line trace marked in a, showing apparent heights of two molecular domains. c, Raw TERS spectra extracted from a single pixel in the datacube marked in the right STM image simultaneously acquired during TERS imaging (–0.1 V, 1 nA, 77 nm2, 3232 pixels, 1 s/pixel). The boundary between the molecular domains is highlighted with a yellow line.d, Reconstructed TERS images based on single-peak analysis for the Raman peak at ~701 cm–1, 1000 cm–1, 1220 cm–1, 819 cm–1, 918 cm–1, 998 cm–1. Theintegration ranges are indicated below the images. The intensity contrasts in these images are too low to make a conclusive chemical identification based on the analysis of weak or overlapped spectral peaks.

S3. Data pre-processing methods for multivariate analysis

The TERS imaging datacube was first pre-processed with denoising, baseline correction, and vector normalization. Denoising was accomplished using a total-variation constraint with a Lagrange parameter of 7 between the data fidelity term and the total variation term,running for 40 iterations. The baseline was modeled as a 7th-order polynomial, and corrected by Asymmetric Least Squares (AsLS) with an asymmetry parameter p=0.05. Vector normalization means that the integral of each spectrum is set to be equal to 1.

Supplementary Figure 3 | TERS imaging datasets in PCA sub-spaces and corresponding dendrograms and merge distance evaluation graphs for molecular domains (a, b and c) and molecular chains (d, e and f), respectively.

S4. Brief description of vertex component analysis(VCA)

As discussed above and in the main text, the HCA has a drawback, in which each spectrum must be assigned to one cluster or another. However, the nanocavity plasmon field in the tunneling gap shows some spatial broadening in real space, which violates the assumption underlying HCA, especially for the TERS spectra acquired at the boundary of two different adjacent molecules.

Supplementary Figure 4| Schematic drawing of VCA.a, Spectral mixtures of 3 pure components form a simplex in a 2-dimensional space. b,Dependence of simplex dimensionality on the number of mixing components. c, TERS imaging dataset of molecular domains in a PCA sub-space, showing mixing between two components along the line defined by the PC1 axis. Purest pixels (endmembers) are at the extrema of the data clouds, marked with small red circles. d, The VCA endmember spectra for two endmembers identified in c. e, The VCA 2-endmember image, with the abundances at points P1, P2, and P3 given in the table below the image.

Vertex component analysis (VCA) overcomes the limitation in HCA by assuming that each measurement is a linear combination of contributions from several underlying components (called “endmembers”). The basic idea of VCA is that if each spectrum in dataset is a mixture of N pure components, all points in the dataset should lie in the interior of an (N-1)-dimensional simplex, where the vertices of the simplex are the pure components and the interior points are varying mixtures of these components. Thus, this simplex should have the smallest volume possible while still enclosing all the points.

As an example, FigureS4a shows a dataset made of 3 components. In the absence of noise, all possible mixtures of the 3 components must lie in the shaded triangle (the mixing simplex). Depending on how many components are in the mixture, the mixing simplex will have different dimensionality, as shown in FigureS4b.Taking the TERS imaging on molecular domains as an example, at the first step we presumethat there are two mixing components in the datacubeconsidering that there are two species in the measured area, therefore the mixing simplex is a line (the PC1 axis), as shown in FigureS4c. In the lower left, a color-coded histogram of the dataset onto the PC 1 axis is also shown, in which the “purest” pixels lie at the extrema of the PC1 axis,and mixed pixels lie towards the center. The “purest” pixels are the VCA endmemberswhose spectra are shown in FigureS4d. The dataset can then be fitted to these endmembers with non-negative least squares fitting (NNLS) to determine the endmember abundance (contribution of each endmember to the pixel spectrum). The resulted abundance image is shown in FigureS4e. Three example pixels are selected to show the endmember abundance in each pixel, showing the mixing.

HCA and VCA were chosen for data evaluation in our analysis due to their simplicity and the correspondence between their assumptions and the high-resolution TERS mapping experimental data. Compared to other clustering methods such as k-means, HCA is deterministic and fully maps the dataset’s connectivity via the dendrogram. The dendrogram is a key component of our analytical pipeline, as its shape can help to determine the likely number of spectral components. Of other modern methods of spectral unmixing and spectral estimation, the most common competitor to VCA is multivariate curve resolution (MCR), particularly the alternating least squares implementation (MCR-ALS). While MCR-ALS has several advantages compared to other methods, including imposing non-negativity constraints on both concentrations and recovered spectra, its implementation is complex and its performance in real situations is often similar to or worse than VCA.1 VCA is extremely simple to implement, requiring the user to only input the dataset and estimated endmember number. Its main constraint is that the dataset should include at least one “pure” pixel for each endmember. Due to the unprecedented resolution of TERS, this requirement is easily satisfied in the datasets described in our experiments.

S5. Estimation ofthe spatial resolution ofTERS imaging

Supplementary Figure 5| Estimation of the spatial resolution of TERS imaging on molecular domains.a, the 2-endmember VCA abundance image. b and c are the point spread functions (PSFs) estimation based on two endmembers abundance profile indicated by a green line in a. Black curves show raw abundance data of ZnTPP (b) and H2TBPP (c) endmembers from green line in a. Dashed blue line and dashed magenta line show the denoised abundance profiles, which approximate edge spread functions. The solid blue and magenta lines indicate the derivatives of these curves, which approximate PSFs. The red curve shows a Gaussian fit to each PSF, with the FWHM indicating an approximateresolution of ~0.4nm.

The spatial resolution is a very important parameter for TERS experiments and is limited by both the pixilation of the images as determined by the step size, as well as the spatial broadening of the localized enhancement of the electric field that produces the TERS spectra. The VCA abundance profiles have substantially higher SNR than the raw images produced by the single-peak analysis shown earlier in FigureS2. In the dataset of ZnTPP and H2TBPP molecular domains, the boundary between the two domains constitutes a step-edge that can be used to explore the spatial resolution of the system. This is similar to a “knife-edge” test common in standard optical imaging. As shown in FigureS5, a line profile across the domain boundary shows a step edge.The derivative of the line profile approximates a point spread function (PSF). A common shorthand for the resolution of an optical system is the full width at half maximum (FWHM)of this PSF. Based on a Gaussian fit of our experimental PSF, we can approximate the resolution of our system as being about 0.4 nm. Note that the pixilation of our system is 0.22 nm/pixel in this experiment, meaning that based on the Nyquist criterion the pixel-limited resolution would be twice this pixilation value, namely 0.44nm. Thus, the resolution here is actually pixel limited, and could potentially be improved even beyond this value by recording data at a smaller lateral step size.

S6. Simulation of TERS-VCAimages for a ZnTPP domain on Ag(111)

The theoretical simulations of the TERS spectra are mainly based on our previous model proposed for single porphyrin molecules (Figure S6a).2, 3 According to Ref. 2, the molecular configuration of a ZnTPP adsorbed on Ag(111) is sketched in Figure S6b, in which the phenyl substituents show a twist angle of roughly 1° and a bending angle of roughly 5° relative to the porphyrin plane. The site-specific TERS spectra of a ZnTPP molecule are simulated by the density functional theory (DFT) calculations via the Gaussian 09 software with the B3LYP/6-31G(d) basis4. To account for the local response within a single molecule under the plasmonic tip, the tip is treated as a highly localized electrical field that is assumed to follow a Gaussian beam distribution with a FWHM of ~0.41 nm (a spatial resolution determined in Supplementary Section S5). The details of this calculation method can be found in our previous publications.2, 3

Supplementary Figure 6 | Simulation of TERS-VCA images for a ZnTPP domain on Ag(111). a, ZnTPP configuration defined in terms of two angles (α, β). Bending angle α quantifies the out-of-plane bending of the phenyl group and twist angle β describes the rotation of the phenyl group about the connecting σ-bond to the porphyrin core. The molecular orientation is defined by the three Euler angles (φ, θ, ψ) according to the molecular coordinates (x, y, z) relative to the STM coordinate (X, Y, Z), with the tilting angle θ being most critical.b, Top view (upper panel) and side view (lower panel) of the ZnTPP configuration adsorbed on Ag(111) used in theoretical simulations.The configuration parameters are listed below the lower panel. c, Theoretical simulation of the TERS-VCAabundance image for an individual unit cell of the ZnTPP domain, with the lattice constants given on the right.dand e,Simulated VCA abundance image (d) and associated endmember spectra (thick lines) for a ZnTPP molecular domain, with experimental endmember spectra (thin lines) also plotted for comparison (e).

In order to simulate the TERS spectra on a ZnTPP domain,a periodic boundary condition (m=1.5 nm, n=1.4 nm, and θ=94°) is introduced to take into account the contributions from the eight nearest neighbor molecules to the TERS signals, as shown in Figure S6c and S6d. Specifically, for a certain position r of the local field, the total TERS intensity can be expressed as , in which denotes the TERS intensity of a single molecule for the local field position r, and Nnear=1 for the nearest neighboring molecules. By scanning the local field over the unit cell (32×32 pixels), a series of simulated TERS spectra can be obtained.

Then, a similar VCA procedure to that describedin Supplementary Section S4 was applied to these simulated TERS spectra of the unit cell by setting the number of endmembers to two (Figure S6c).By shifting the unit cell with lattice constant (m=1.5 nm, n=1.4 nm, and θ=94°), the TERS-VCAabundance image on a ZnTPP domain is obtained (Figure S6d).The associated endmember spectra are also shown in Figure S6e.

We would like to note that the simulation procedure used in the present work is based on a simplified model without considering molecule–substrate and intermolecular interactions. The model used in our simulation is too simple to give a good agreement with the experimental spectra. With such simple model, however, the periodic pattern in the TERS-VCA image on the ZnTPP domain is fairly well reproduced. The use of such a simple model is mainly to understand qualitatively the physical origin for the periodic pattern in the TERS-VCA image on the ZnTPP domain. More sophisticated calculations are needed to improve the agreement with experimental spectra.

S7. Single-peak analysisof TERS images on a molecular chain

Supplementary Figure 7 | Single-peak analysis of TERS images on a molecular chain. a,TERS spectra extracted from the datacube that are averaged over the blue, green and grey squares (3×3 pixels) marked in the right STM image simultaneously acquired during TERS imaging (–0.3 V, 0.8 nA, 1616 nm2, 3232 pixels, 1 s/pixel).b, Reconstructed TERS images based on single-peak analysis.The peak positions and theassociated integration ranges are indicated above and below the images, respectively. Again, the intensity contrasts in these images are too low to make a conclusive chemical identification for two different molecules, based on the analysis of weak or overlapped spectral peaks.

S8. Hierarchical clustering analysis (HCA) of TERS images on a molecular chain

Supplementary Figure 8 | HCA images on molecular chains with different cluster numbers and corresponding cluster average spectra.

For the VCA in TERS imaging on the molecular chains, a direct application of VCA to the dataset leads to an incorrect estimate of the background spectrum. This is because the featureless background can be considered as the absence of any spectrum rather than an endmember of its own, and thus violates the assumptions of VCA. To overcome this, we first use 2-cluster HCA to split the dataset into background and molecular regions, and exclude the background regions from the VCA analysis. Once the VCA endmembers of the molecular regions are determined, the spectral fitting uses those VCA-determined endmembers plus the cluster average spectrum of the background as the pure components in the NNLS fitting process.

S9. References

1.Zhang X, Tauler R. Application of Multivariate Curve Resolution Alternating Least Squares (MCR-ALS) to remote sensing hyperspectral imaging. Anal Chim Acta 2013; 762: 25-38.

2.Zhang R, Zhang Y, Dong ZC, Jiang S, Zhang C et al. Chemical mapping of a single molecule by plasmon-enhanced Raman scattering. Nature 2013; 498: 82–86.

3.Jiang S, Zhang Y, Zhang R, Hu C, Liao M et al. Distinguishing adjacent molecules on a surface using plasmon-enhanced Raman scattering. Nat Nanotechnol 2015; 10: 865–869.

4.Frisch MJ, Trucks GW, Schlegel HB, Scuseria GE, Robb MA et al. Gaussian 09. Wallingford, CT, USA: Gaussian, Inc.; 2009.