# Determining by Raman spectroscopy the average thickness and N-layer-specific surface coverages of MoS<sub>2</sub> thin films with domains much smaller than the laser spot size

Felipe Wasem Klein, Jean-Roch Huntzinger, Vincent Astié, Damien Vairy, Romain Parret, Houssin Makhlouf, Sandrine Juillaguet, Jean-Manuel Decams, Sylvie Contreras, Périne Landois, et al.

► To cite this version:

Felipe Wasem Klein, Jean-Roch Huntzinger, Vincent Astié, Damien Vairy, Romain Parret, et al. Determining by Raman spectroscopy the average thickness and N-layer-specific surface coverages of MoS<sub>2</sub> thin films with domains much smaller than the laser spot size. Beilstein Journal of Nanotechnology, 2024, 15, pp.279-296. 10.3762/bjnano.15.26. hal-04747945

HAL Id: hal-04747945

https://hal.science/hal-04747945v1

Submitted on 28 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Determining by Raman spectroscopy the average thickness and *N*-layer-specific surface coverages of MoS₂ thin films with domains much smaller than the laser spot size

Felipe Wasem Klein¹, Jean-Roch Huntzinger¹, Vincent Astié², Damien Voiry³, Romain Parret⁴, Houssine Makhlouf¹, Sandrine Juillaguet¹, Jean-Manuel Decams², Sylvie Contreras¹, Périne Landois¹, Ahmed-Azmi Zahab¹, Jean-Louis Sauvajol¹ and Matthieu Paillet*¹

Full Research Paper

Open Access

Address:

¹Laboratoire Charles Coulomb, Université de Montpellier, CNRS, F-34095, Montpellier, France, ²Annealsys, 139 Rue des Walkyries, 34000 Montpellier, France, ³Institut Européen des Membranes, IEM, UMR 5635, Université Montpellier, ENSCM, CNRS, Montpellier, France and ⁴Aix Marseille Université, CNRS, CINAM, UMR 7325, Campus de Luminy, 13288, Marseille, France

Email:

Matthieu Paillet* - matthieu.paillet@umontpellier.fr

Keywords:

molybdenum disulfide; number of layers; Raman spectroscopy; thin film; transition metal dichalcogenides

Beilstein J. Nanotechnol. 2024, 15, 279-296.

https://doi.org/10.3762/bjnano.15.26

Received: 13 November 2023

Accepted: 20 February 2024

Published: 07 March 2024

This article is part of the thematic issue "Applications of Raman spectroscopy in the characterization of nanomaterials".

Guest Editor: P. Araujo

<figure><img src="image_3.png" /></figure>

© 2024 Wasem Klein et al.; licensee

Beilstein-Institut.

License and terms: see end of document.

## Abstract

Raman spectroscopy is a widely used technique to characterize nanomaterials because of its convenience, non-destructiveness, and sensitivity to materials change. The primary purpose of this work is to determine via Raman spectroscopy the average thickness of MoS₂ thin films synthesized by direct liquid injection pulsed-pressure chemical vapor deposition (DLI-PP-CVD). Such samples are constituted of nanoflakes (with a lateral size of typically 50 nm, i.e., well below the laser spot size), with possibly a distribution of thicknesses and twist angles between stacked layers. As an essential preliminary, we first reassess the applicability of different Raman criteria to determine the thicknesses (or layer number, *N*) of MoS₂ flakes from measurements performed on reference samples, namely well-characterized mechanically exfoliated or standard chemical vapor deposition MoS₂ large flakes deposited on 90 ± 6 nm SiO₂ on Si substrates. Then, we discuss the applicability of the same criteria for significantly different DLI-PP-CVD MoS₂ samples with average thicknesses ranging from sub-monolayer up to three layers. Finally, an original procedure based on the measurement of the intensity of the layer breathing modes is proposed to evaluate the surface coverage for each *N* (i.e., the ratio between the surface covered by exactly *N* layers and the total surface) in DLI-PP-CVD MoS₂ samples.

<footnote>* Corresponding author</footnote>
Introduction

The advent of two-dimensional (2D) layered materials beyond
graphene has initiated a new field of research [1-3]. In the
family of 2D layered structures, transition metal dichalco-
genides (TMDs) have attracted considerable attention from
academia and regarding potential applications [4-9] because of
a number of remarkable properties [10-12]. Particularly, it was
found that the properties of layered TMDs drastically change
when their thickness is reduced to a monolayer [13,14]. Lay-
ered TMD structures have a graphite-like structure with each
graphene sheet replaced with an X–M–X or MX₂ triatomic
layer, where X is a chalcogen atom (e.g., sulfur, selenium, or
tellurium) and M is a transition metal atom (e.g., molybdenum
or tungsten) [10].

Among the layered TMD materials, molybdenum disulfide, MoS₂, is of particular interest in optoelectronic applications because of its transition to a direct bandgap semiconductor with very high photoluminescence quantum yield when thinned down to a monolayer [13-17]. Its unique electronic and optical properties could provide an edge in many future applications.

The multilayers MoS₂ structures are of the most common 2Hc
type, where atomic layers are arranged in such way that the
stacking between two adjacent layers corresponds to a twist
angle of θ = 60°, and any Mo atom is sitting on top of two S
atoms of the adjacent layers [18,19]. However, during the syn-
thesis process (e.g., chemical vapor deposition (CVD) synthe-
sis) or when using precise transfer or AFM tip manipulation
techniques [20], twisted MoS₂ can be formed with two adjacent
layers stacked with a relative twist angle (θ) varying from 0 to
60°. Such twisted-layered MoS₂ structures can exhibit a variety
of interesting physical properties including unconventional
super conductivity [21,22], non-linear optics [23,24], and moiré
excitons [25].

Because the properties of MoS₂ flakes are first a function of
their thickness, or layer number (N), it is of a primary impor-
tance to determine the *N* of MoS₂ flakes, including twisted
MoS₂ flakes and defective MoS₂ flakes, synthesized by differ-
ent ways. Independently of the structural organization between
adjacent layers, a MoS₂ flake is usually named *NL-MoS₂*, or
simply *NL*, with *N* being the number of MoS₂ triatomic layers,
which defines the thickness of the flake.

Several optical techniques have been developed to identify the
*N* of MoS₂ flakes produced by different methods. Among these
techniques, Raman spectroscopy is widely used thanks to its
convenience, non-destructiveness, and sensitivity to materials
change, including strain, temperature, doping, and defects [26].
Concerning the characterization of MoS₂ flakes, different infor-

mation can be derived from the measurement of the Raman fea-
tures (frequencies, linewidths, and intensities) of intralayer
phonon modes as well as those of the interlayer modes, the
so-called layer breathing (LB) modes and shear (S) modes.

Recently, we have developed the reproducible direct growth of wafer-scale MoS₂ thin films on SiO₂/Si substrates by direct liquid injection pulsed-pressure chemical vapor deposition (DLI-PP-CVD) using low-toxicity precursors [27]. Such MoS₂ thin films showed good stoichiometry (Mo/S = 1.94–1.95) and the potential for high photoluminescence quantum yield. However, atomic force microscopy revealed that they are constituted of nanoflakes (with a lateral size of typically 50 nm) with possibly a distribution of thicknesses. Furthermore, depending on the synthesis conditions, the MoS₂ surface coverage can be incomplete, and the thin film average thickness can vary. These samples thus have characteristics, especially thickness inhomogeneities smaller than the laser spot size, that differ from the ones used to establish Raman spectroscopy-based MoS₂ layer counting methods [26,28-33]. In this context, the primary purpose of this work is to develop and validate an approach for determining the average thickness of such sub-laser spot size inhomogeneous MoS₂ thin films using Raman spectroscopy.

First, we reassess here as a ground work the information that
can be derived from the Raman spectra of MoS₂ flakes for the
evaluation of their thickness, *N*. Different Raman criteria for the
determination of the thicknesses of MoS₂ flakes are first recol-
lected; after the specification of the experimental protocol,
domains and limits of application of these criteria are precisely
defined from measurements performed on reference samples.
These samples are well-characterized, either mechanically
exfoliated or standard CVD MoS₂ large flakes deposited on
90 ± 6 nm SiO₂ on Si substrates. Then, we determine which
Raman information is relevant to estimate the average thick-
ness of MoS₂ samples produced by the DLI-PP-CVD method,
which are constituted of nanoflakes and, thus, significantly dif-
ferent from the reference samples. Finally, an original proce-
dure based on the layer breathing mode intensities is proposed
to evaluate the surface coverage for each *N*, that is, the ratio be-
tween the surface covered by exactly *N* layers and the total sur-
face, in DLI-PP-CVD samples.

Results and Discussion
Experimental procedure

To define a robust experimental Raman protocol to evaluate the
thickness of a MoS₂ flake (i.e., its number of layers, *N*), it is
first necessary to specify some parameters that can have a direct
influence on the quality and accuracy of the results. The first
parameter is the wavelength of the incident laser light used in
the Raman experiments. As it will be detailed in the following, the measurements of the frequencies, linewidths, and intensity of first-order Raman active phonon modes of MoS₂ have to be obtained with good accuracy in order to evaluate the thickness of a MoS₂ flake. These phonons modes are (i) the in-plane phonon mode involving relative motion of Mo and S atoms with E' symmetry for a monolayer (E¹₂g for bulk) and (ii) the out-of-plane phonon mode involving only out-of-plane motions of S atoms with A'₁ symmetry for a monolayer (A₁g for bulk). These modes are located around 385 and 405 cm⁻¹, respectively, in neutral and defect-free MoS₂ monolayers [33,34]. More precisely, in MoS₂ multilayers, the symmetries of these phonon modes are E' and A'₁ for an odd number of layers, and Eg and A₁g for an even number of layers. For simplicity, hereafter when we will discuss the dependence on *N* of the features of these phonon modes, they will be simply referred to as in the bulk, E¹₂g and A₁g, independently of the number of layers.

A drastic change of the Raman spectra, especially in the fre-
quency range of the A₁g and E¹₂g modes, occurs when the spec-
tra are excited at an energy close to those of the A and B exci-
tons located around 655 nm (1.89 eV) and 601 nm (2.06 eV),
respectively, in MoS₂ monolayers [35,36]. When the incident
laser energy is in the range of the A and B exciton energies (the
so-called resonance conditions), other bands associated to dif-
ferent second-order processes are observed in the Raman spec-
tra with a strong intensity, their frequencies, widths, and intensi-
ty depending on the excitation energy [36]. In addition, reso-
nance conditions alter the symmetry selection rules of phonons
of MoS₂ [35]. Some of the second-order bands overlap with the
A₁g and E¹₂g modes, complicating the exact determination of
the parameters of these modes recorded under resonance condi-
tions. Furthermore, since the MoS₂ exciton characteristics
(energy, width, and spectral weight) can be changed by several
factors (e.g., stacking, strain, doping, and defects), the Raman
intensities measured with a single laser wavelength close to
exciton energies can be affected by external factors and differ
for samples elaborated by different methods. For these reasons
and in the aim to use Raman spectroscopy to count the number
of MoS₂ layers, one must necessarily work under off-resonance
conditions, that is, by using incident laser energy far from both
exciton resonance energies. In this work, we chose to perform
Raman experiments using 532 nm (2.33 eV) laser excitation,
because this is sufficiently far from the energy range of A and B
excitons [35].

All Raman spectra reported in this paper were recorded on dif-
ferent samples deposited on SiO₂/Si(100) substrate. Hence, the
second parameter essential to define is the SiO₂ thickness.
Indeed, the multiple reflection interferences that occur in the
air/MoS₂/SiO₂/Si structure influence significantly the intensity

of the phonon modes [28,37]. In this work, we chose to focus
on substrates with a SiO₂ thickness around 90 nm, which corre-
sponds to the first optimum value for MoS₂ monolayer (N = 1)
Raman enhancement with a 532 nm excitation energy and also
amplifies the signal in the wavelength range of photolumines-
cence emission (around 650 nm).

The third parameter to define is the power of the 532 nm light, Pλ, impinging the sample. Much of the Raman information available to evaluate the thickness of MoS2 flakes is based on the following parameters: (i) on precise measurements of frequency of the A1g and E12g phonon modes of MoS2. These lead to a precise knowledge of the frequency difference ΔωA-E. It was established that ΔωA-E depends monotonously on the number of layers, and ΔωA-E is largely used as criterion to evaluate the thickness of MoS2 flakes [26,29,30]; (ii) on the precise evaluation of the integrated intensities of the phonon modes of MoS2, namely A(A1g) and A(E12g), with respect to the integrated intensity of the 521 cm⁻¹ mode from a bare area of the oxidized silicon substrate, A₀(Si), used as an intensity reference [31], or from the silicon substrate underneath the MoS₂ flake, A₂D(Si) [28]; (iii) on the precise measurement of the A₂D(Si)/A₀(Si) intensity ratio [31]; and (iv) on the measurement of ultralow-frequency modes, the so-called breathing modes and shear modes. The frequencies and the number of LB and S modes allow one to identify the number of layers [32,33] and the presence of twist between adjacent layers from the vanishing of the S modes in twisted MoS₂ flakes [20,38-41].

Then, it is essential to determine the limit value of the laser
power so that the above measurements are not affected by laser
irradiation. Figure 1 shows the evaluation of the temperature of
MoS2 flakes prepared in different ways and that of the Si sub-
strate as functions of the laser power impinging on the sample
through a 100× objective (N.A. 0.9). The power was cycled
between ≈5 µW and ≈2 mW. The temperature of MoS2 flakes
is evaluated from the Stokes/anti-Stokes intensity ratio of A1g
phonon modes (similar results are obtained using E12g) and that
of silicon from the Stokes/anti-Stokes intensity ratio of the
521 cm⁻¹ Si mode (see [42] for method details). While the
silicon temperature is quasi-insensitive to Pλ, the temperature
of MoS2 flakes changes monotonically, reversibly, and quasi-
linearly with Pλ (see inset of Figure 1). For MoS2, we found an
increase rate of 25–30 °C/mW for monolayers (1L-MoS2) and
40–45 °C/mW for bilayers (2L-MoS2). Usual effects of sample
heating are the frequency shift of the phonon modes and their
concomitant broadening. In Supporting Information File 1, the
frequency and the linewidth of the Si mode as functions of the
laser power are displayed (Figure S2). These two parameters are
found to be insensitive to Pλ below 0.5 mW. More intriguing is
the evolution of the frequency (Figure 2a) and width (Figure 2b)
<figure><img src="image_4.png" /><figcaption>**Figure 1:** Evolution of the temperature of Si substrate (black circles) and MoS₂ flakes (1L exfoliated: green upward triangles, 1L CVD: blue downward triangles, and 2L CVD: red squares) as functions of the incident laser power in log-scale. The absolute values of the temperature under the laser spot were measured from the variation of the Stokes/AntiStokes ratio of the 521 cm⁻¹ Si mode for the substrate and of the A₁g phonon modes for MoS₂ flakes. Inset, same data plotted with a linear laser power scale.</figcaption></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 2:** CVD 1L-MoS₂. Evolution of A'₁ (red squares) and E' (blue circles) Raman modes frequencies (a) and full width at half maximum (b) as functions of the incident laser power during a cycle from 5 µW up to 2 mW and back to 5 µW. (c) Evolution of the A'₁ relative shift versus the E' relative shift as a function of laser power. The color code of each point corresponds to the incident laser power as displayed on the color bar. The data are compared to the expected evolutions for heating effect only [43] (red line) and for doping only [44] (magenta line). To enable direct comparison, the strain contribution has been removed, and the corresponding zero doping point is labeled as well as the directions corresponding to *p* and *n* doping.</figcaption></figure>

of the phonon modes of 1L-MoS₂ as a function of $P_{\lambda}$. For a simple thermal effect [43] and given the 25–30 °C/mW temperature increase rate determined previously, the frequencies of A'₁ and E' modes should both downshift by 0.3–0.4 cm⁻¹/mW and the width of A'₁ should increase by ≈0.2 cm⁻¹/mW (the width of E' should remain constant) contrary to what is observed in Figure 2a,b. To clarify this point, we present in Figure 2c (filled dots) the relative shift of the frequency of the A'₁ mode versus that of the E' mode measured on CVD 1L-MoS₂ at different laser powers. In the same plot the expected shifts of these modes are reported (i) as functions of a pure thermal effect (Figure 2c, red line) [43] and (ii) as functions of the doping state (Figure 2c, magenta curve) [44]. Clearly the relative shift of the A'₁ mode frequency versus that of the E' mode frequency as a function of $P_{\lambda}$ significantly differs from the behavior expected by considering a simple thermal effect. Consequently,

the results reported in Figure 2c clearly evidence photo-doping of 1L-MoS₂ concomitant with a thermal effect, as already observed for MoS₂ on SiO₂/Si [45] as well as for graphene [42]. Furthermore, the evolution of the A'₁ and E' widths with $P_{\lambda}$ (Figure 2b), that is, the weak change of the E' width and the significant increase of the A'₁ width concomitant with the A'₁ frequency decrease, support this interpretation [44]. For a laser power smaller than 0.3 mW, photo-doping remains rather low, but it is the dominant contribution to the shift of the modes. Similar results were obtained on other samples including exfoliated 1L-MoS₂. The effects were found irreversible in some cases when $P_{\lambda}$ exceeded 1 mW but remained always reversible if $P_{\lambda}$ was kept below 1 mW.

Based on the above information, all Raman results reported and discussed in this paper were obtained by using a $P_{\lambda}$ around
0.1 mW chosen as a good compromise between mitigating laser effects and maintaining measurement efficiency to ensure the accuracy of the Raman criteria discussed in the next part of this paper.

In summary, unless specified otherwise, all Raman spectra reported and discussed in this paper were recorded at an excitation wavelength of 532 nm, with a laser working-power close to 0.1 mW, and using a 100× objective (N.A. 0.9), on MoS₂ flakes or thin films deposited on SiO₂/Si substrate with a SiO₂ thickness of 90 ± 6 nm.

## Application of Raman criteria to characterize MoS₂ flakes

In this part, we report and discuss the advantages and limits of some Raman criteria that were found to be efficient to derive the thickness (i.e., the number of layers *N*) of large MoS₂ flakes prepared by different ways, namely mechanical exfoliation and standard CVD (including twisted CVD 2L-MoS₂). Then, we discuss the application of Raman spectroscopy to characterize samples synthesized by DLI-PP-CVD. In contrast to the first two kinds of MoS₂ samples, the latter are constituted of nanoflakes with possibly a distribution of thicknesses and twist angles between adjacent layers of multilayer domains as well as a higher number of defects.

## Exfoliated MoS₂ flakes as reference samples

We performed Raman experiments on mechanically exfoliated MoS₂ [1] that will serve as reference samples. The stacking sequence in exfoliated MoS₂ flakes is of the 2Hc-type [34]. The common feature of all these samples is to have a limited number of defects. Note also that all exfoliated flakes have a lateral size (few micrometers at minimum) significantly larger than the diameter of the laser spot. In such flakes, the exact number of layers, *N*, is determined by combining optical microscopy, spectral reflectivity, and the measurement of the breathing modes and shear modes in the ultralow frequency (ULF) range of the spectra [32-34].

One of the most popular criteria to determine the number of layers of MoS₂ flakes is the measurement of $\Delta\omega_{A-E}$, that is, the frequency difference between the $A_{1g}$ and $E_{2g}$ phonon modes [26,29,30]. Figure 3a shows the dependence of $\Delta\omega_{A-E}$ on the number of layers measured on exfoliated MoS₂ flakes deposited on Si/SiO₂ substrates with four different SiO₂ thicknesses. As previously well documented in the literature, we confirm that $\Delta\omega_{A-E}$ depends monotonously on the number of layers and does not depend on the SiO₂ thickness (Figure 3a). The separation between *N* and *N* + 1 values are larger than the experimental uncertainties (error bars in the graph) up to *N* = 3. The error bars start to overlap between *N* = 4 and *N* = 5. Comparison with

data from the literature (see inset in Figure 3a) shows that this overlap occurs even between *N* = 3 and *N* = 4 when additional variability due to setup and samples is taken into account. Above *N* = 4, the separation becomes too small compared to the uncertainty. Thus, the measurement of $\Delta\omega_{A-E}$ in exfoliated MoS₂ flakes allows one to evaluate with good accuracy the number of layers for *N* ≤ 3. It is then necessary to supplement the $\Delta\omega_{A-E}$ criterion with others to reliably count thick multilayers. In addition, we will establish in the following that this criterion has to be taken with care to derive *N* in MoS₂ samples other than reference exfoliated MoS₂, because the $A_{1g}$ and $E_{2g}$ frequencies, and thus the value of $\Delta\omega_{A-E}$, can be affected by different factors such as stacking order, strain, doping, and defects which can be present in MoS₂ flakes prepared by other ways [44,46-49].

To evaluate the number of layers, we can also use information associated with the integrated intensity of MoS₂ phonon modes. Figure 3b and Figure 3d show, respectively, the dependences of the normalized integrated intensities of the $A_{1g}$ and the $E_{2g}$ mode as functions of *N* for four values of the SiO₂ thickness. For normalization, we use here an external reference, which is a bare Si(111) wafer with only native oxide. In the following, $A(\text{Si}_{111})$ stands for the integrated intensity of the Si(111) 521 cm⁻¹ mode. This reference is preferred to the Si(100) substrate with 90 ± 6 nm SiO₂ to avoid the effects of the SiO₂ thickness variations and crystal orientation. For comparison with other setups or references, the polarization ratio of our setup and the relative values measured on Si(100) with native oxide and 90 nm SiO₂ are given in Supporting Information File 1. As noted by several authors and predicted by the optical interference model, the normalized integrated intensities of MoS₂ modes, namely $A(\text{A}_{1g})/A(\text{Si}_{111})$ and $A(\text{E}_{2g})/A(\text{Si}_{111})$, increase first with *N* and then decrease showing a maximum for *N* = 4–5 for all SiO₂ thicknesses. Obviously, this non-monotonous dependence prevents using these measurements alone to evaluate the number of layers for *N* > 4. Moreover, a significant dependence of the MoS₂ Raman intensity on the SiO₂ thickness occurs for *N* > 2, pointing out the importance to determine precisely this latter parameter.

Another criterion to derive the thickness of MoS₂ flakes is the $A_{2D}(\text{Si})/A_0(\text{Si})$ intensity ratio [31]. For the evaluation of this ratio, it is of great practical advantage to use the same silicon (the silicon below the oxide, which is Si(100) in the present work) in the measurement of $A_{2D}(\text{Si})$ and $A_0(\text{Si})$. A necessary precaution is that the Si(100) substrate orientation has to be kept the same for both measurements. Another advantage is to give a common origin to the plots of $A_{2D}(\text{Si})/A_0(\text{Si})$ as a function of $N$ ($A_{2D}(\text{Si})/A_0(\text{Si}) = 1$ for $N = 0$) for any SiO₂ thickness.
<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 3:** Mechanically exfoliated samples on SiO₂/Si substrates with four different SiO₂ thicknesses specified in the legends. Dependence on the number of MoS₂ layers (a) of the frequency difference between the A₁g and E¹₂g phonon modes, and of the normalized integrated intensities (see text) of (b) the A₁g, (c) the Si 521 cm⁻¹, and (d) the E¹₂g modes. The inset in (a) shows a comparison of the average values measured here with data from the literature. The insets in (b–d) show the corresponding average values for each SiO₂ thickness and each number of layers.</figcaption></figure>

Figure 3c displays the A₂D(Si)/A₀(Si) ratio measured on exfoliated MoS₂ flakes deposited on SiO₂/Si(100) substrates with four different SiO₂ thicknesses as a function of N. We confirm the monotonous decrease of the A₂D(Si)/A₀(Si) ratio with increasing N [31]. For each N, the value of this ratio depends on the SiO₂ thickness (Figure 3c; black, blue, green, and red symbols correspond to a SiO₂ thickness of 84, 87, 89, and 96 nm, respectively). Despite the monotonous dependence of the A₂D(Si)/A₀(Si) ratio, its dependence on SiO₂ thickness combined with experimental errors lead to the conclusion that the measured values for N and N + 1 can overlap for any N if the SiO₂ thickness is not known with good accuracy. For a given SiO₂ thickness, the gap between the A₂D(Si)/A₀(Si) ratio for N and N + 1 is sufficient to ensure a rather good reliability only for N ≤ 5.

In summary, for exfoliated MoS₂, considering jointly the three Raman criteria (i) value of ΔωA-E, (ii) value of the normalized integrated intensities of the A₁g and E¹₂g modes, and (iii) value of the A₂D(Si)/A₀(Si) ratio, one can unambiguously derive the number of layers as long as N ≤ 4 and the SiO₂ thickness is precisely known. It has also been suggested in the literature

to use the intensity ratio A(A₁g)/A₂D(Si) (or equivalently A(E¹₂g)/A₂D(Si)). As it will be discussed in the following, we see two major problems with this approach. The first is the dependence of the Si signal on the crystal orientation and the SiO₂ thickness. The second relates to the fact that using this ratio, rather than using each measurand independently and contrasting them, even if more practical, can hide some information.

Finally, we compare in Supporting Information File 1, Figure S3, the dependence on N of A(A₁g)/A(Si₁₁₁) and A₂D(Si)/A₀(Si) for three SiO₂ thicknesses and two microscope objectives with different numerical apertures, N.A. = 0.9 (blue symbols in Figure S3), and N.A. = 0.5 (red symbols in Figure S3). We observe that the normalized integrated intensity of the A₁g mode, A(A₁g)/A(Si₁₁₁) (Figure S3a–c), is independent of the value of N.A. Concerning the dependence on N of A₂D(Si)/A₀(Si) (Figure S3d–f), we found that this ratio is slightly smaller for N.A. = 0.5 than for N.A. = 0.9 and is in a better agreement with the model of Li and coworkers [31] (black solid line in Figure S3d–f). (Note that in this latter work the experimental data on which the model has
been adjusted were recorded using a numerical aperture N.A. ≈ 0.45).

## MoS₂ flakes prepared by CVD

In this part, we analyze the pertinence of the previous criteria to derive the thickness of large MoS₂ flakes synthesized by CVD. In a first part, we probe the effectiveness of these criteria to evaluate the thickness of large standard CVD MoS₂ flakes. Such flakes have a limited number of defects and, like in exfoliated MoS₂, the stacking sequence is of the 2Hc type. In the second part, we examine the relevance of these criteria to evaluate the thickness of twisted CVD MoS₂ flakes.

**Standard CVD MoS₂ flakes:** As derived from the features of the LB and S ultralow-frequency modes, these samples do not show any twist between adjacent layers (presence of S modes for all flakes with $N \ge 2$). The flakes are thus characterized by 2Hc stacking (or close to 2Hc stacking) and a low number of defects, and are named standard CVD MoS₂ flakes. On the

basis of the latter features, the structure of these flakes is close to that of exfoliated MoS₂ flakes. However, the high temperature used in the CVD synthesis and interaction with the substrate can lead to lattice distortion and the presence of vacancies and doping. In the following, we limit our study to a number of layers $N \le 4$.

Figure 4a compares the values of $\Delta\omega_{A-E}$ measured on exfoliated (Figure 4a, black symbols) and standard CVD MoS₂ flakes (Figure 4a, red symbols) for $N \le 4$. As previously, the exact number of layers is obtained by combining optical microscopy, spectral reflectivity, and number and frequencies of LB and S modes. For both kinds of MoS₂ flakes, $\Delta\omega_{A-E}$ increases monotonously with $N$, but for the same $N$, the values of $\Delta\omega_{A-E}$ are systematically larger in standard CVD MoS₂ flakes. We attribute this discrepancy mainly to a difference of strain states between the two kinds of samples. Exfoliated samples are mostly found with low or slight compressive strain, while CVD samples are under tension. Other factors such as doping, defects,

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

**Figure 4:** Comparison between CVD (red filled dots) and mechanically exfoliated (black-line-connected open squares) MoS₂ on a 87 nm SiO₂/Si substrate. Dependence on the number of MoS₂ layers of (a) the frequency difference between the $A_{1g}$ and $E^{1}_{2g}$ phonon modes, and of the normalized integrated intensities (see text) of (b) the $A_{1g}$, (c) the Si 521 cm⁻¹, and (d) the $E^{1}_{2g}$ modes. For the 4L CVD, only few points are measured and the errors are not shown because they cannot be properly derived.
or stacking were shown to lead to large changes of $\Delta\omega_{A-E}$ [44,46,49]. In summary, the value of $\Delta\omega_{A-E}$ is clearly and significantly sample-dependent. Consequently, $\Delta\omega_{A-E}$ cannot be considered as a definitive criterion to derive the number of layers in any of MoS$_2$ flakes prepared in different ways. In other words, one cannot define a single master curve, $\Delta\omega_{A-E}$ vs N, which would be valid for all the MoS$_2$ flakes independently of their preparation method or environment.

The dependencies on N of the normalized integrated intensities of $A_{1g}$ and $E_{2g}$ modes and the $A_{2D}(Si)/A_0(Si)$ ratio measured on standard CVD flakes are compared with the average values of exfoliated samples with the same substrate $SiO_2$ thickness (Figure 4b–d). In contrast to $\Delta\omega_{A-E}$, the N dependencies of these intensities are very close in exfoliated and standard CVD MoS$_2$ flakes. Only $A_{2D}(Si)/A_0(Si)$ and the normalized integrated intensity of the $E_{2g}$ phonon modes slightly differ for

$N=4$. However, this may be due to the fact that the statistics is rather poor on this measurement, because this flake is rather small compared to the others. With regards to these results, these measurands give important information to evaluate the number of layers of 2Hc-stacked MoS$_2$ flakes independently of the elaboration procedure as long as $N \le 4$.

**Twisted CVD MoS$_2$ flakes:** Other interesting samples are large CVD MoS$_2$ flakes that present a twist angle, $\theta$, between adjacent layers. We exemplify here the complexity to characterize such samples from the previous Raman criteria with the case of twisted MoS$_2$ bilayers. The identification of the bilayer character of the investigated flakes was unambiguously obtained independently from spectral reflectivity and optical contrast.

Figure 5a shows the low-frequency range of spectra recorded on three types of MoS$_2$ bilayer (named 2L-MoS$_2$ in the following),

<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

**Figure 5:** $\theta$-2L-MoS$_2$. (a) Raman spectra of 3R-, 2Hc- and 30°-2L-MoS$_2$ (black, blue, and red lines, respectively) in the ultralow-frequency region. Layer breathing (LB) and shear (S) modes are labelled. (b) Raman spectra of 2Hc-, 20°- and 30°-2L-MoS$_2$, 1L-MoS$_2$, and bare substrate in the region of 130–240 cm$^{-1}$. Raman modes with $\theta$-dependent frequencies (FTA and FLA) are labelled. (c) Raman spectra of 2Hc-, 3R-, 20°- and 30°-2L-MoS$_2$, and 1L-MoS$_2$ in the region of 360–430 cm$^{-1}$. $E_{2g}$, $A_{1g}$, and $\theta$-dependent $FA'_1$ are labelled. (d) Raman spectra of 2Hc-, 3R-, 20°- and 30°-2L-MoS$_2$, 1L-MoS$_2$, and bare substrate in the region of the Si 521 cm$^{-1}$ mode. Y-scale intensities are normalized to the incident laser power and the acquisition time. Presented spectra are averages extracted from Raman maps.
namely a bilayer with $\theta \approx 30^\circ$ (this sample belongs to the so-called twisted-bilayer family for which $0 < \theta < 60^\circ$ and is named in the following as θ-2L-MoS₂), the so-called 2Hc-2L-MoS₂ and the so-called 3R-2L-MoS₂. In the latter structure, the stacking between two adjacent layers corresponds to a twist angle of $\theta = 0^\circ$, and it is such that the S atoms of the top monolayer are superimposed on the Mo atoms of the bottom monolayer, and the Mo atoms of the top monolayer are above the hexagon centers of the bottom monolayer [50]. The spectrum in the low-frequency range is dominated by the contributions of the LB and S modes, the frequencies of these modes depending on the twist angle [20,39]. The LB mode emerges in the Raman spectra of all 2L-MoS₂ samples. In line with previous results [20], the Raman shift of the peak position of the LB mode in 3R-2L-MoS₂ is smaller than that of the 2Hc-2L-MoS₂, and the LB mode Raman shift in 30°-2L-MoS₂ is even smaller. Also, in agreement with the literature [20,38-41], the S mode vanishes in 30°-2L-MoS₂.

In Figure 5b, the 130–240 cm⁻¹ range of the Raman spectra recorded on monolayer, 2Hc-2L-MoS₂, and two θ-2L-MoS₂ is displayed. As well documented in the literature, this frequency range is dominated by the contributions of second-order Raman processes [20,26]. The general profile of the spectra is similar in 1L-MoS₂, 2Hc-2L-MoS₂, and θ-2L-MoS₂ with the exception that in the latter flakes, new bands, named FLA and FTA, are superimposed to the second-order Raman spectra. The FLA and FTA modes in θ-2L-MoS₂ are attributed, respectively, to folded longitudinal acoustic phonons and folded transverse acoustic phonons of the monolayer due to the presence of a moiré superlattice [20]. As shown in the literature [20], the frequencies of these modes depend on the twist angle (see Figure 2e in [20]). Unfortunately, these dependencies show a mirror behavior with respect to $\theta = 30^\circ$. This means that from given FLA and FTA positions, two values are possible: $\theta \in [0,30]^\circ$ or its mirror $60^\circ - \theta$. As a consequence, $\theta$ will be given in the range of 0–30° in all plots in Figure 6 with the possibility that the values attributed to θ-2L-MoS₂ could be $60^\circ - \theta$ instead. For instance, the data from both 2Hc-2L-MoS₂ and 3R-2L-MoS₂ are reported at $\theta = 0^\circ$ in these plots. From the positions of the FTA and FLA, we claim that the spectra of the two θ-2L-MoS₂ displayed in Figure 5 correspond to 20°-2L-MoS₂ (Figure 5b–d, solid green line) and 30°-2L-MoS₂ (Figure 5a–d, solid red line), respectively.

The dependence of $A_{1g}$ and $E^{1}_{2g}$ modes on the twist angle (derived from the positions of FTA and FLA modes) is reported in Figure 5c. The frequency of the $E^{1}_{2g}$ mode in 3R-, 2Hc-, and θ-2L-MoS₂ is downshifted with respect to its frequency in 1L-MoS₂, and it does not show a clear dependence on the twist angle. In contrast, the profile of the $A_{1g}$ mode signifi-

cantly depends on the twist angle, and a new mode, named FA'₁, appears on the high-frequency side of the $A_{1g}$ mode. The FA'₁ mode is identified as Raman scattering from moiré phonons associated with the A'₁ dispersion curve of 1L-MoS₂. It is folded onto the zone center and, consequently, becomes Raman active [20]. Obviously, its frequency depends on the twist angle and the θ-dependence of the FA'₁ frequency was recently established both theoretically and experimentally (see Figure 3e in [20]). On the basis of these previous results, we have been able to evaluate the value of θ for each 2L-MoS₂ investigated from the position of the FA'₁ mode. The values of the angles derived from the position of FTA/FLA and FA'₁ are in close agreement.

The objective of this work is to characterize the thickness of all MoS₂ flakes. The relevance of the criteria based on the frequency (Δω<sub>A-E</sub>) and normalized integrated intensity (A(A<sub>1g</sub>)/A(Si<sub>111</sub>)) of the A<sub>1g</sub> mode has to be reevaluated in twisted 2L-MoS₂ flakes. As shown in Figure 5d, the normalized intensity of the 521 cm<sup>-1</sup> Si mode from the substrate underneath MoS₂ flakes is close in all the 2L-MoS₂ and independent of the twist angle.

Figure 6 summarizes and details the dependence on the twist angle of the four Raman criteria defined above for 2L-MoS₂. In all plots of Figure 6, the values of angles were derived from the positions of FTA, FLA and FA'₁. The values of the different criteria measured for θ-2L-MoS₂ are compared with the average values of the same criteria measured on exfoliated 1L-, 2L-, and 3L-MoS₂ flakes. In θ-2L-MoS₂, the value of Δω<sub>A-E</sub> strongly depends on the twist angle and significantly differs from the average value measured in 2Hc-2L-MoS₂ (Figure 6a). For θ = 30°, the value of Δω<sub>A-E</sub> is close to the one found in CVD 1L-MoS₂ [51]. In consequence, using Δω<sub>A-E</sub> alone could lead to a wrong evaluation of the thickness of twisted 2L-MoS₂.

The normalized integrated intensity A(A₁g)/A(Si₁₁₁) significantly decreases when the twist angle increases (Figure 6b, red dots), and in 30°-2L-MoS₂, the value of A(A₁g)/A(Si₁₁₁) is close to the average value found in 1L-MoS₂ (Figure 6b, blue solid line). The behavior of A(A₁g)/A(Si₁₁₁) is opposite to the one of the normalized integrated intensity of the FA'₁ mode, A(FA'₁)/A(Si₁₁₁), the latter increasing with the twist angle (Figure 6b, gray squares). These results are in qualitative agreement with those reported in [40]. It can be emphasized that the integrated intensity of A₁g and FA'₁ bands taken together (Figure 6b, orange triangles) is close to the average value found for 2Hc-2L-MoS₂ (Figure 6b, red solid line). The reason for this compensation between A(A₁g) and A(FA'₁) is not clear yet, but it could present a practical advantage in the use of the global integrated intensity of the spectral band, located around the po-
<figure><img src="image_20.png" /><figcaption>**Figure 6:** θ-2L-MoS₂. θ-dependence of the frequency difference between the A₁g and E¹₂g phonon modes (a); θ-dependence of the normalized integrated intensities (see text) of A₁g and FA'₁ modes and their sum (b), of the Si 521 cm⁻¹ mode (c), and of the E¹₂g mode (d). In (b–d), the corresponding average values for 1L-MoS₂, 2Hc 2L-MoS₂, and 3L-MoS₂ are plotted as horizontal lines (blue, red, and green, respectively) for comparison.</figcaption></figure>

sition of the A₁g mode for the evaluation of the thickness of twisted MoS₂ flakes.

We also observed a tendency for A(E¹₂g)/A(Si₁₁₁) to be slightly lower for θ-2L-MoS₂ than for 2Hc-2L-MoS₂ (or 3R-2L-MoS₂, which is similar), but to a lesser extent compared to A(A₁g)/A(Si₁₁₁), that is ca. 20% vs ca. 50% at maximum, respectively (Figure 6d). These results are also in qualitative agreement with those reported in [40]. Finally, only the value of the A₂D(Si)/A₀(Si) ratio seems to provide a robust/reliable information to characterize the thickness of MoS₂ flakes, since it is found largely independent of θ in all measured 2L-MoS₂ samples (Figure 6c). Even if further work is needed to complete the data presented here with other values of θ and twisted MoS₂ samples with N > 2, we anticipate that the value of A₂D(Si)/A₀(Si) ratio would be close in twisted and 2Hc-stacked MoS₂ multilayers. However, as previously recalled, the sensitivity of this ratio to the SiO₂ thickness and the gap between the A₂D(Si)/A₀(Si) ratios for N and N + 1 permit to ensure the determination of N with a rather good reliability only for N ≤ 5.

## DLI-PP-CVD MoS₂ nanoflakes

The aim of this part is to define which Raman information is relevant to estimate the thickness of MoS₂ samples produced by DLI-PP-CVD. These samples are significantly different from the previous ones (exfoliated and standard CVD). Indeed, they are constituted of nanoflakes (with a lateral size of typically 50 nm, i.e., well below the laser spot size) with possibly a distribution of thicknesses and twist angles between adjacent layers of multilayer domains and a higher number of defects (the average inter-defect distance ranges from 3 to 6 nm as estimated from the LA and A₁g intensity ratio [52]). In addition, the MoS₂ surface coverage is a priori unknown and can be incomplete. It is then necessary to implement a first check criterion that ensures that the thickness estimation method based on the comparison with results obtained on exfoliated samples is still valid. More generally, this point is critical for the characterization of samples synthesized using new methods or new precursors that can lead to the co-deposition of several by-products (such as carbon, oxides, and metals), which can significantly change the measured Raman intensities. Based on the results presented in the previous sections, we have shown that the value of the
$A_{2D}(\text{Si})/A_0(\text{Si})$ ratio provides a robust/reliable Raman information to characterize the thickness of $\text{MoS}_2$ flakes for $N \le 5$. However, this parameter does not rely unambiguously on the presence of $\text{MoS}_2$. The deposition of any other material would influence its value and could lead to a wrong estimation. In the most general case, the sample characteristics are not perfectly known and can be significantly different from the reference characteristics. As a consequence, it seems mandatory to compare the thickness estimated from the $A_{2D}(\text{Si})/A_0(\text{Si})$ ratio with other measurands directly related to the presence of $\text{MoS}_2$. To this aim, we propose to use jointly the normalized integrated intensity of the $\text{MoS}_2$ phonon modes, namely $A(A_{1g})/A(\text{Si}_{111})$ and/or $A(E^{1}_{2g})/A(\text{Si}_{111})$.

In Figure 7, the values of $A(A_{1g})/A(\text{Si}_{111})$ (Figure 7a) and $A(E^{1}_{2g})/A(\text{Si}_{111})$ (Figure 7b) are plotted as functions of the value of $A_{2D}(\text{Si})/A_0(\text{Si})$. In these graphs, the data obtained on DLI-PP-CVD samples are compared with the average reference measurements established previously on exfoliated $\text{MoS}_2$ deposited on Si/SiO$_2$ substrate with the same SiO$_2$ thicknesses, namely 96 nm (red open dots in Figure 7) and 87 nm (blue open dots in Figure 7). Note that in the exfoliated samples, the exact number of layers $N$ is perfectly known and given on the plots of Figure 7 close to corresponding open dots. The idea behind this representation comes from the expectation that the presence of contaminations or deposition of other species would have a different impact on the Raman intensity coming from $\text{MoS}_2$ in the film and on the one coming from the Si substrate underneath the deposited thin film. It is, thus, expected that the measurements on contaminated or highly defective $\text{MoS}_2$ thin films will fall off the reference curve. Indeed, data obtained on poorly

crystalline $\text{MoS}_2$ films synthesized by DLI atomic layer deposition (not shown) are found systematically and significantly below the corresponding reference curve. Concerning the DLI-PP-CVD samples presented in Figure 7, the $A(A_{1g})/A(\text{Si}_{111})$ vs $A_{2D}(\text{Si})/A_0(\text{Si})$ dependence is found fully compatible with the respective reference exfoliated curves (Figure 7a). The $A(E^{1}_{2g})/A(\text{Si}_{111})$ vs $A_{2D}(\text{Si})/A_0(\text{Si})$ data points mainly agree for thin layers ($A_{2D}(\text{Si})/A_0(\text{Si}) > 0.8$) but fall systematically below the corresponding reference exfoliated curves for thicker layers ($0.8 > A_{2D}(\text{Si})/A_0(\text{Si}) > 0.6$) as shown in Figure 7b.

Another way to compare the results is estimating the thickness of DLI-PP-CVD samples by interpolation from exfoliated data of the measured values for $A_{2D}(\text{Si})/A_0(\text{Si})$, $A(A_{1g})/A(\text{Si}_{111})$, and $A(E^{1}_{2g})/A(\text{Si}_{111})$. In Figure 8a (respectively 8b), the average number of layers ($\bar{N}$) obtained using $A(A_{1g})/A(\text{Si}_{111})$ (respectively $A(E^{1}_{2g})/A(\text{Si}_{111})$) are plotted as a function of the number derived from $A_{2D}(\text{Si})/A_0(\text{Si})$. It can be emphasized that non-integer values are found for $\bar{N}$, indicating the presence of a mix with unknown proportions of bare substrate (0L), 1L-$\text{MoS}_2$, 2L-$\text{MoS}_2$, 3L-$\text{MoS}_2$, and so on in the investigated DLI-PP-CVD films. It is also noticeable that the errors of $\bar{N}$ estimated from $A(A_{1g})/A(\text{Si}_{111})$ become larger when $\bar{N}$ is close to 3 as a consequence of the smoother dependence of this parameter with $\bar{N}$. In agreement with the conclusion drawn above from Figure 7, Figure 8a illustrates the coherence between the values of $\bar{N}$ derived from $A(A_{1g})/A(\text{Si}_{111})$ and $A_{2D}(\text{Si})/A_0(\text{Si})$. All data remain close to the red solid line that represents the ideal relation $y[\bar{N}$ via $A(A_{1g})/A(\text{Si}_{111})] = x[\bar{N}$ via $A_{2D}(\text{Si})/A_0(\text{Si})]$. Figure 8b as well confirms that the values of $\bar{N}$ derived from $A(E^{1}_{2g})/A(\text{Si}_{111})$ and $A_{2D}(\text{Si})/A_0(\text{Si})$ agree well for $\bar{N} < 1.5$, but

<figure><img src="image_21.png" /><figcaption>**Figure 7:** DLI-PP-CVD samples. Normalized integrated intensities of (a) the $A_{1g}$ and (b) the $E^{1}_{2g}$ modes as functions of the normalized integrated intensities of the Si 521 cm$^{-1}$ mode of DLI-PP-CVD samples on 87 nm (green squares) and 96 nm (black squares) SiO$_2$ on Si substrates. Each point corresponds to the average value extracted from 121-point Raman maps. The average values measured on corresponding exfoliated $\text{MoS}_2$ samples are also plotted for comparison as open dots (blue and red, respectively, for 87 nm and 96 nm SiO$_2$), and the corresponding number of layers are provided. The lines are guides to the eye.</figcaption></figure>
<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /><figcaption>**Figure 8:** DLI-PP-CVD samples. Average number of layers of DLI-PP-CVD samples obtained by interpolation from exfoliated data of normalized integrated intensities of (a) the A₁g and (b) the E¹₂g modes, and of (c) microreflectivity spectra as functions of the average number of layers obtained by interpolation from exfoliated data of normalized integrated intensities of the Si 521 cm⁻¹ mode. The red lines in (a–c) correspond to y = x. (d) Frequency difference between the A₁g and E¹₂g phonon modes of DLI-PP-CVD samples as a function of the average number of layers obtained by interpolation from exfoliated data of normalized integrated intensities of the Si 521 cm⁻¹ mode. Green (respectively black) squares correspond to DLI-PP-CVD samples synthesized on 87 nm (respectively 96 nm) SiO₂/Si substrates. In (d), the average values measured on corresponding exfoliated (respectively CVD) MoS₂ samples are also plotted for comparison as open blue (respectively red) dots. The lines are guides to the eye.</figcaption></figure>

the values of $\bar{N}$ from $A(E^{1}_{2g})/A(\text{Si}_{111})$ are systematically lower than those obtained from $A_{2D}(\text{Si})/A_{0}(\text{Si})$ when $\bar{N} > 1.5$. One explanation could the presence of a larger proportion of multilayer regions in the thicker samples, for which, if they are twisted, $A(E^{1}_{2g})/A(\text{Si}_{111})$ has been shown to be attenuated in the previous section. If so, the question then arises why the same behavior is not observed for $A(\text{A}_{1g})/A(\text{Si}_{111})$ contrary to what would be expected. A possibility could be that because of the observed broadening of the A₁g mode in DLI-PP-CVD samples (presumably due to local heterogeneities in terms of doping, strain, defects, or thickness), the FA'₁ mode becomes indistinguishable from the A₁g mode. As a consequence, the intensity of FA'₁ would merge with A(\text{A}_{1g}) and compensate its attenuation. Other explanations relying on the presence of defects or strain cannot be disregarded, and further works are needed to fully clarify this point.

In order to further confirm the validity of the estimations of $\bar{N}$ for DLI-PP-CVD samples, we compare in Figure 8c the $\bar{N}$

values derived from $A_{2D}(\text{Si})/A_{0}(\text{Si})$ with the ones obtained independently from spectral microreflectivity. A very good agreement is found between the two series of data, establishing definitively the relevance of the $A_{2D}(\text{Si})/A_{0}(\text{Si})$ ratio to give with good accuracy the average thickness of DLI-PP-CVD MoS₂ samples for $\bar{N} \le 3$. This agreement justifies the use of the values of $\bar{N}$ derived from $A_{2D}(\text{Si})/A_{0}(\text{Si})$ as abscissa axis in the previous plots.

Finally, in Figure 8d the frequency difference between the A₁g and E¹₂g phonons is plotted as a function of $\bar{N}$ estimated from $A_{2D}(\text{Si})/A_{0}(\text{Si})$ for DLI-PP-CVD samples and compared to the data obtained on exfoliated and CVD MoS₂. DLI-PP-CVD data are distributed between the two curves obtained from the reference samples. This further confirms that this measurement cannot be used to evaluate with good accuracy their average thicknesses. Nevertheless this comparison can be informative, showing that samples with $\bar{N} < 1$ are most certainly mainly composed of 1L-MoS₂ and suggesting that the proportions of
2L-MoS₂, 3L-MoS₂, or more gradually increase with $\bar{N}$, which is compatible with AFM observations (not shown).

To get further insight on the number of layer distributions in DLI-PP-CVD samples, we have measured their ULF modes. Representative ULF spectra are shown in Figure 9a for samples with average thicknesses ranging from 0.6 up to 2.8 MoS₂ layers as estimated from A₂D(Si)/A₀(Si). Up to $\bar{N}$ = 1.3, only the LB mode of 2L-MoS₂ is observed around 40 cm⁻¹ [20,38-41], showing that these samples can only be composed of 1L-MoS₂ and twisted 2L-MoS₂ plus possibly uncovered (bare substrate) regions. For thicker samples, the S mode of 2L-MoS₂ around 24 cm⁻¹ is additionally visible, as well as a signal between 25 and 30 cm⁻¹, corresponding to the LB and S modes of 3L-MoS₂ [32,33]. For $\bar{N}$ ≥ 2.4, the LB mode of 4L-MoS₂ is also present around 21 cm⁻¹; there may also be a weak signal around 17 cm⁻¹ (corresponding to the LB mode of 5L-MoS₂) reflecting the presence of 5L-MoS₂. The S mode of 4L-MoS₂ could be present as well around 28 cm⁻¹, but it is hardly distinguishable from the LB and S modes of 3L-MoS₂. Thus, ULF Raman spectra give valuable qualitative information on the different N present in each sample. Quantitative information relies on the determination of the surface coverages for each N ($\sigma_N$), that is, the ratio between the surface covered by exactly N layers and the total surface. With N = 0 standing for the bare substrate and $N_{max}$ being the largest number of layers present in the sample, the definition of the average number of layers $\bar{N}$ can be written as

$$ \bar{N} = \sum_{N=1}^{N_{max}} \sigma_N \times N \quad (1) $$

and the total coverage (including bare substrate areas) is obviously 100%:

$$ \sum_{N=0}^{N_{max}} \sigma_N = 1 \quad (2) $$

AFM imaging (see Supporting Information File 1, Figure S4) reveals that for $\bar{N}$ > 1.25, the surface is fully covered by MoS₂, that is, $\sigma_0$ = 0, which removes an unknown. In addition, for $\bar{N}$ < 1.3, there is no signature of more than two layers, and we can set $\sigma_{N≥3}$ = 0 with confidence. Hence, for 1.25 < $\bar{N}$ < 1.3, the set of Equations 1 and 2 simplifies to

$$ \bar{N} = \sigma_1 \times 1 + \sigma_2 \times 2 \\ \sigma_1 + \sigma_2 = 1 $$

This allows one to readily determine the two remaining unknowns $\sigma_1$ and $\sigma_2$, since $\bar{N}$ is known from A₂D(Si)/A₀(Si).

Hereafter, a linear relationship between the Raman signal $I_{LB}^{2L}$, which is the LB mode peak intensity of 2L (N = 2) areas (the broad but well-identified 40 cm⁻¹ peak), and the surface coverage is assumed, namely $\sigma_2 = \alpha_2 I_{LB}^{2L}$. The ratio $\alpha_2 = \sigma_2/I_{LB}^{2L}$ is determined from five samples (1.25 < $\bar{N}$ < 1.3) for which we now have both the coverage $\sigma_2$ and the Raman signal $I_{LB}^{2L}$.

Because $\alpha_2$ is now known and assuming that the linearity between $\sigma_2$ and the Raman signal $I_{LB}^{2L}$ holds (which should be a good approximation for the thin multilayers considered here), $\sigma_2 = \alpha_2 I_{LB}^{2L}$ can be obtained directly for all samples from the Raman spectra, and is thus no longer an unknown.

<figure><img src="image_26.png" /><figcaption>**Figure 9:** (a) Raman spectra in the ultralow-frequency region of selected DLI-PP-CVD samples with different average thicknesses as labelled in the legend. Presented spectra are averages extracted from 121-point Raman maps. (b) Surface coverages of bare substrate (black open squares), 1L-MoS₂ (red open dots), 2L-MoS₂ (blue open up-triangles), 3L-MoS₂ (green open down-triangles), and 4L-MoS₂ (orange open diamonds) as functions of the average DLI-PP-CVD sample thickness. The model results are shown as full lines. See the main text and Supporting Information File 1 for the description of the model.</figcaption></figure>
We now turn to the samples with $\bar{N} < 1.25$, which may present some bare substrate areas, so $\sigma_0$ and $\sigma_1$ are a priori unknown. Both $\bar{N}$ and $\sigma_2$ are determined as explained above, and $\sigma_{N\ge 3} = 0$ is again a safe estimate. Hence, Equations 1 and 2 reduce to the system

$$ \begin{aligned} \bar{N} &= \sigma_1 \times 1 + \sigma_2 \times 2 \\ \sigma_0 + \sigma_1 + \sigma_2 &= 1 \end{aligned} $$

which may be solved trivially for $\sigma_0$ and $\sigma_1$.

A similar approach can be used for samples with $1.3 < \bar{N} < 2$ that are fully covered ($\sigma_0 = 0$) and might present some trilayers but show no trace of thicker layers. We set $\sigma_{N\ge 4} = 0$, and the system of equations reduces to

$$ \begin{aligned} \bar{N} &= \sigma_1 \times 1 + \sigma_2 \times 2 + \sigma_3 \times 3 \\ \sigma_1 + \sigma_2 + \sigma_3 &= 1 \end{aligned} $$

Both $\bar{N}$ and $\sigma_2$ are determined as explained above. Thus, the system can be trivially solved for the two remaining unknowns $\sigma_1$ and $\sigma_3$.

At this point it would be natural to get the proportionality between $\sigma_3$ and a Raman signal attributed to 3L areas ($N = 3$), and proceed recursively to obtain $\sigma_4$ in slightly thicker layers, and so on and so forth. In practice this becomes challenging because of the uncertainty on the 3L ($N = 3$) Raman signal, which is less clear than the 2L ($N = 2$) peak. Another approach gave better results.

Three samples with $\bar{N}$ between 2.75 and 2.85 are thick enough to neglect $\sigma_0$ and $\sigma_1$, yet thin enough for $\sigma_5$ to also be negligible as a first approximation. Equations 1 and 2 then reduce to

$$ \begin{aligned} \bar{N} &= \sigma_2 \times 2 + \sigma_3 \times 3 + \sigma_4 \times 4 \\ \sigma_2 + \sigma_3 + \sigma_4 &= 1 \end{aligned} $$

where $\bar{N}$ and $\sigma_2$ are known, so $\sigma_3$ and $\sigma_4$ can be determined readily.

The LB mode of 4L-MoS$_2$, located around $21 \text{ cm}^{-1}$, is sufficiently separated from other modes to be identified (which was not the case for $I_{\text{LB}}^{3L}$), so that $I_{\text{LB}}^{4L}$ can be extracted from the spectra. From the three $2.75 < \bar{N} < 2.85$ samples $\alpha_4 = \sigma_4/I_{\text{LB}}^{4L}$ is determined. Then, assuming again a linear relationship $\sigma_4 = \alpha_4 I_{\text{LB}}^{4L}$, the coverage by 4L ($N = 4$) layers can be determined for all samples. This removes another unknown.

Now the last remaining case of $2 < \bar{N} < 2.75$ samples can be solved, as

$$ \begin{aligned} \bar{N} &= \sigma_1 \times 1 + \sigma_2 \times 2 + \sigma_3 \times 3 + \sigma_4 \times 4 \\ \sigma_1 + \sigma_2 + \sigma_3 + \sigma_4 &= 1 \end{aligned} $$

give $\sigma_1$ and $\sigma_3$ directly, since $\bar{N}$, $\sigma_2$ and $\sigma_4$ are known.

The results obtained using this procedure are shown in Figure 9b where $\sigma_N$ (with $N$ from 0 to 4) is plotted as a function of $\bar{N}$, the average sample thickness. On this graph, all values of $\sigma_2$ (respectively $\sigma_4$) are calculated using $I_{\text{LB}}^{2L}$ (respectively $I_{\text{LB}}^{4L}$) even for the samples used to derive the proportionality coefficient $\alpha_2$ (respectively $\alpha_4$). As shown in Figure 9b for samples with $1.25 < \bar{N} < 1.3$ (respectively $2.75 < \bar{N} < 2.85$), we find by this way $-0.03 < \sigma_0 < 0.04$ (respectively $-0.01 < \sigma_1 < 0.01$) with little fluctuations around the expected value of 0.

Just below the full coverage of the sample surface by MoS$_2$ ($\sigma_0 > 0$), both $\sigma_1$ and $\sigma_2$ increase with a slight tendency of $\sigma_2$ to increase faster than $\sigma_1$. Indeed, 1L-MoS$_2$ represents 80-90% of the deposited MoS$_2$ for $\bar{N} = 0.5$ and 70-80% for $\bar{N} = 1.3$. The maximum of $\sigma_1$ is reached around $\bar{N} = 1.3$ when the sample surface is totally covered by MoS$_2$ ($\sigma_0 = 0$), and $\sigma_1$ starts to decrease above this value. Around $\bar{N} = 1.6$, 1L-MoS$_2$ only represents 50% of the MoS$_2$. $\sigma_2$ continues to increase and reaches a maximum value of ≈50% around $\bar{N} = 2$ and then decreases for thicker samples. 3L-MoS$_2$ starts to appear after the substrate surface is completely covered by MoS$_2$ and increases continuously, representing about 50% of the thickest samples ($\bar{N} \approx 2.8$).

In order to verify our approach, we implemented a 2D growth toy model (see Supporting Information File 1 for details). The model results are shown in Figure 9b as full lines and give a good agreement with the experimental results. It should be noted that within this representation ($\sigma_N = f(\bar{N})$), the results of the model are remarkably robust to any parameter changes (the curves are almost insensitive to either doubling or halving the cell size and, thus, the advance rate, or to multiplying or dividing the growth rate by 5). In other words, this means that this comparison with the experiment cannot be used to validate any model parameters but demonstrates the relevance of the proposed procedure to estimate the $\sigma_N$ from the experiments. Nevertheless, it has to be noted that while for $\bar{N} < 1.3$ the ULF Raman signature of 2L-MoS$_2$ remains very similar, it is not the case for thicker samples with the notable appearance of the S mode of 2L-MoS$_2$ around $24 \text{ cm}^{-1}$ [32,33]. This could mean that the stacking order distribution changes. As a consequence, the hypothesis based on the proportionality between $I_{\text{LB}}^{2L}$ and $\sigma_2$ would probably be less valid above $\bar{N} = 1.3$, and an error on the absolute values deduced can be expected. However, the appearance of the S mode of 2L-MoS$_2$ around $24 \text{ cm}^{-1}$ could also be related to NL-MoS$_2$ (with $N \ge 3$) constituted of a stacking sequence.
quence where 2L are not twisted, for example, the so-called
t(1+2)L, t(2+2)L, ... structures [26,53]. In this case, our hypoth-
esis would remain more appropriate. Despite this unknown as
well as the other approximations made, we believe that the main
tendencies can still be captured by the proposed analysis.
Further works are needed to determine and improve the accu-
racy of the method.

Conclusion

In this work we have reviewed all Raman information leading
to the evaluation of the thickness of MoS₂ flakes, that is, the
layer number *N*. First, we have analyzed in detail the effects of
some experimental parameters, namely the wavelength of the
incident laser light used in the experiments, the power of the
incident light, and the oxide thickness of the SiO₂/Si substrate
on which the flakes are deposited, on the quality and accuracy
of Raman results. Based on this analysis, an experimental
protocol has been defined and systematically applied to large
MoS₂ flakes (i.e., single-domain flakes much larger than the
laser spot), including twisted MoS₂ flakes, prepared by differ-
ent methods on the one hand and to MoS₂ thin films composed
of nanoflakes prepared by the DLI-PP-CVD method on the
other hand. Special attention was paid to the measurement
statistics.

The limits of different Raman criteria which allow one to deter-
mine the thicknesses of MoS2 flakes, namely (i) the value of
Δω<sub>A-E</sub>, (ii) the value of the normalized integrated intensity
of A<sub>1g</sub> and E<sup>1</sup><sub>2g</sub> MoS<sub>2</sub> modes, and (iii) the value of the
A<sub>2D</sub>(Si)/A<sub>0</sub>(Si) ratio, have been precisely studied in the differ-
ent types of MoS<sub>2</sub> samples. We definitely confirm that Δω<sub>A-E</sub>
cannot be considered a robust criterion to derive the number of
layers in MoS<sub>2</sub> samples. We found that the value of the
A<sub>2D</sub>(Si)/A<sub>0</sub>(Si) ratio provides the most robust/reliable informa-
tion to characterize the thickness of MoS<sub>2</sub> large flakes, espe-
cially since it is found largely independent of the twist angle.
The limit of application of this criterion is N ≤ 5, under the
condition that the SiO<sub>2</sub> thickness is precisely known.

We then apply this analysis procedure to DLI-PP-CVD sam-
ples constituted of nanoflakes with a lateral size of typically
50 nm (well below the laser spot size) with possibly a distribu-
tion of thicknesses and twist angles between adjacent layers of
multilayer domains and a higher number of defects. Our results
definitively establish the relevance of the A₂D(Si)/A₀(Si) ratio to
give with good accuracy their average thickness N̄, for N̄ ≤ 3.
Nevertheless, we emphasize that this criterion is not only
related to the presence of MoS₂ and can be influenced by
several factors, such as the co-deposition of by-products or the
presence of defects, leading to a wrong estimation of N̄. We
propose to combine A₂D(Si)/A₀(Si) with the normalized inte-

grated intensity of the MoS₂ phonon modes, namely A(A₁g) and/or A(E¹₂g). Although limiting the application to N̄ ≤ 3, this approach enables the validation of the A₂D(Si)/A₀(Si) ratio to determine N̄ in the presented case, and we anticipate that it would avoid possible errors in unfavorable situations.

Finally, to get further insight on the number of layer distribu-
tions in DLI-PP-CVD samples, we have measured their ULF
modes. An original procedure based of the measurement of
the intensity of the layer breathing modes allows one to eval-
uate the surface coverage (σ<sub>N</sub>) for each N. A 2D growth toy
model gives a good agreement with the experimental results
supporting the proposed procedure to estimate the σ<sub>N</sub> from the
ULF spectra.

Experimental

Samples preparation

Mechanical exfoliation

MoS₂ flakes were obtained by micromechanical cleavage of a MoS₂ crystal (HQ graphene) using scotch tape (Nitto) and PDMS slabs (Gel-pak). They were then transferred onto Si substrates with SiO₂ layers of different thicknesses, namely 84, 87, 90, and 96 nm. Flakes were selected by optical microscopy and their thicknesses were determined by optical contrast.

Standard CVD process

MoS₂ was grown by CVD on 87 nm SiO₂ on Si substrates using MoO₃ (Sigma-Aldrich, 25 mg) and sulfur (Sigma-Aldrich, 250 mg) powders as solid precursors using a 1 inch quartz tube furnace. MoO₃ powder was placed in the center of the heating zone of the furnace, while sulfur was placed upstream at the furnace inlet. Prior to growth, air was evacuated by flowing Ar (ultrahigh purity, Linde) for 15 min at 200 sccm, after which the tube was heated to 200 °C for 10 min. The temperature was then increased to 750 °C under Ar (100 sccm), and it was held at this value for 15 min before cooling naturally to room temperature.

Direct-liquid injection pulsed-pressure chemical vapor deposition (DLI-PP-CVD)

The 12 × 11 mm SiO₂/Si (with 87 nm or 96 nm SiO₂ thick-
nesses) substrates were cleaned in acetone (C₃H₆O, technical,
Acros Organics), isopropanol (C₃H₇OH, 99.8%, Höfer Chemie
GmbH) and deionized water (H₂O, Acros Organics) under
ultrasonic agitation for 10 min each, before being blown dry
with nitrogen. They were then immediately loaded on the
susceptor of the reaction chamber (Annealsys MC-050) for
deposition. Solutions of 0.001 M molybdenum hexacarbonyl
(Mo(CO)₆, 98%, Strem Chemicals) and 0.002 M sulfur (S,
99.999%, Acros Organics) in anhydrous toluene (C₆H₅CH₃,
99.8%, Sigma-Aldrich) were used as precursors. The process is
as follows: Following sample installation, the chamber is closed and brought to about 0.02 mbar. For monolayer depositions, it is imperative that the substrate is thoroughly cleaned and free of adsorbates. Therefore, to ensure complete desorption of remaining contaminants, the samples were kept for 30 min under vacuum at room temperature inside the deposition chamber. For the first part of the process, the pumping direction is reversed so that all species are pumped from the deposition chamber to the back of the reactor.

Nitrogen (800 sccm) is flowed through the chamber (200 sccm through the gas line, and 600 sccm through two injection heads) and the substrates are brought to 750 °C at a ramp of 2 °C/s. The reactor is kept in this state for 5 min for homogenization purposes. While still in reverse direction pumping, 0.3 g/min of both precursors are injected and vaporized to prepare the evaporation system for deposition. Then, the Mo(CO)₆ injection is stopped, the pumping direction is switched back to the deposition direction and hydrogen (40 sccm) is added to the gas mix. For 1 min, sulfur is injected to clean any remaining contaminants, and to prepare the surface of the substrate for MoS₂ deposition, then the deposition works in 20 s cycles. During one cycle, a single pulse of 3 to 10 ms of Mo(CO)₆ is injected while the S injection is set to 0.3 g/min. This 20 s cycle is repeated 80 to 160 times. The quantity of MoS₂ deposited is controlled by the quantity of Mo(CO)₆ injected, that is, the pulse duration and the number of cycles.

## Raman spectroscopy

Raman spectra and maps were recorded using an Acton spectrometer fitted with a Pylon CCD detector and a 1800 grooves/mm grating (≈0.6 cm⁻¹ between each CCD pixel). The samples were excited with a 532 nm (2.33 eV) laser (Newport Millennium Prime or Cobolt Samba) through an Olympus microscope objective either 100× (numerical aperture 0.9) or 50× (numerical aperture 0.5). The full width at half-maximum (FWHM) of the focused laser spot with the 100× objective is about 380 nm. Optimized focus conditions were checked for each measurement. The samples were mounted on a three-axis piezoelectric stage (Physik Instrumente) to ensure the precise positioning and focusing of the laser spot. A Si(111) wafer with only native oxide sample was used as a daily reference for the system. The laser power was continuously measured during acquisitions allowing for intensity normalizations of the Raman spectra at each point of the maps. All data presented in this paper, unless specified otherwise, are extracted from Raman maps constituted by hundreds to thousands points (see Supporting Information File 1 for an example), which were analyzed using a custom-made software. All reported points are the average values obtained by Gaussian fitting of the data distribution extracted from Raman maps (corresponding to

hundreds to several thousands of spectra), and the error bars correspond to 99.7% confidence intervals (±3 standard deviations).

## 2D growth toy model

The model used the DynamicGrids.jl package, which was part of the Dispersal.jl framework [54], see Supporting Information File 1 for more details.

## Supporting Information

Supporting Information File 1 contains additional figures with an example of Raman maps, the Si mode as a function of the laser power, a comparison between two microscope objectives, other intensity references, atomic force microscopy images, and details of the 2D growth toy model. Supporting Information File 2 is a recording of the growth simulation.

**Supporting Information File 1**

Additional experimental data.
[https://www.beilstein-journals.org/bjnano/content/supplementary/2190-4286-15-26-S1.pdf]

**Supporting Information File 2**

Recording of the growth simulation.
[https://www.beilstein-journals.org/bjnano/content/supplementary/2190-4286-15-26-S2.mp4]

## Funding

This work was supported by Occitania region "Recherche et Société(s)" program in the frame of DIMENSION project N°2018-003267.

## ORCID® iDs

Vincent Astié - https://orcid.org/0000-0001-6548-178X
Romain Parret - https://orcid.org/0000-0003-0662-6373
Sylvie Contreras - https://orcid.org/0000-0002-4733-9998
Matthieu Paillet - https://orcid.org/0000-0001-7038-9101

## References

1. Novoselov, K. S.; Jiang, D.; Schedin, F.; Booth, T. J.; Khotkevich, V. V.; Morozov, S. V.; Geim, A. K. *Proc. Natl. Acad. Sci. U. S. A.* **2005**, *102*, 10451–10453. doi:10.1073/pnas.0502848102

2. Gibney, E. *Nature* **2015**, *522*, 274–276. doi:10.1038/522274a

3. Lin, Z.; McCreary, A.; Briggs, N.; Subramanian, S.; Zhang, K.; Sun, Y.; Li, X.; Borys, N. J.; Yuan, H.; Fullerton-Shirey, S. K.; Chernikov, A.; Zhao, H.; McDonnell, S.; Lindenberg, A. M.; Xiao, K.; LeRoy, B. J.; Drndić, M.; Hwang, J. C. M.; Park, J.; Chhowalla, M.; Schaak, R. E.; Javey, A.; Hersam, M. C.; Robinson, J.; Terrones, M. *2D Mater.* **2016**, *3*, 042001. doi:10.1088/2053-1583/3/4/042001
4. Withers, F.; Del Pozo-Zamudio, O.; Mishchenko, A.; Rooney, A. P.; Gholinia, A.; Watanabe, K.; Taniguchi, T.; Haigh, S. J.; Geim, A. K.; Tartakovskii, A. I.; Novoselov, K. S. *Nat. Mater.* **2015**, *14*, 301–306.
doi:10.1038/nmat4205

5. Jariwala, D.; Sangwan, V. K.; Lauhon, L. J.; Marks, T. J.; Hersam, M. C. *ACS Nano* **2014**, *8*, 1102–1120.
doi:10.1021/nn500064s

6. Georgiou, T.; Jalil, R.; Belle, B. D.; Britnell, L.; Gorbachev, R. V.; Morozov, S. V.; Kim, Y.-J.; Gholinia, A.; Haigh, S. J.; Makarovsky, O.; Eaves, L.; Ponomarenko, L. A.; Geim, A. K.; Novoselov, K. S.; Mishchenko, A. *Nat. Nanotechnol.* **2013**, *8*, 100–103.
doi:10.1038/nnano.2012.224

7. Wang, H.; Yu, L.; Lee, Y.-H.; Shi, Y.; Hsu, A.; Chin, M. L.; Li, L.-J.; Dubey, M.; Kong, J.; Palacios, T. *Nano Lett.* **2012**, *12*, 4674–4680.
doi:10.1021/nl302015v

8. Baugher, B. W. H.; Churchill, H. O. H.; Yang, Y.; Jarillo-Herrero, P. *Nat. Nanotechnol.* **2014**, *9*, 262–267. doi:10.1038/nnano.2014.25

9. Zhang, Y. J.; Oka, T.; Suzuki, R.; Ye, J. T.; Iwasa, Y. *Science* **2014**, *344*, 725–728. doi:10.1126/science.1251329

10. Wang, Q. H.; Kalantar-Zadeh, K.; Kis, A.; Coleman, J. N.; Strano, M. S. *Nat. Nanotechnol.* **2012**, *7*, 699–712. doi:10.1038/nnano.2012.193

11. Calman, E. V.; Fogler, M. M.; Butov, L. V.; Hu, S.; Mishchenko, A.; Geim, A. K. *Nat. Commun.* **2018**, *9*, 1895.
doi:10.1038/s41467-018-04293-7

12. Rivera, P.; Schaibley, J. R.; Jones, A. M.; Ross, J. S.; Wu, S.; Aivazian, G.; Klement, P.; Seyler, K.; Clark, G.; Ghimire, N. J.; Yan, J.; Mandrus, D. G.; Yao, W.; Xu, X. *Nat. Commun.* **2015**, *6*, 6242.
doi:10.1038/ncomms7242

13. Mak, K. F.; Lee, C.; Hone, J.; Shan, J.; Heinz, T. F. *Phys. Rev. Lett.* **2010**, *105*, 136805. doi:10.1103/physrevlett.105.136805

14. Splendiani, A.; Sun, L.; Zhang, Y.; Li, T.; Kim, J.; Chim, C.-Y.; Galli, G.; Wang, F. *Nano Lett.* **2010**, *10*, 1271–1275. doi:10.1021/nl903868w

15. Huang, S.; Ling, X.; Liang, L.; Kong, J.; Terrones, H.; Meunier, V.; Dresselhaus, M. S. *Nano Lett.* **2014**, *14*, 5500–5508.
doi:10.1021/nl5014597

16. Scheuschner, N.; Ochedowski, O.; Kaulitz, A.-M.; Gillen, R.; Schleberger, M.; Maultzsch, J. *Phys. Rev. B* **2014**, *89*, 125406.
doi:10.1103/physrevb.89.125406

17. Eda, G.; Yamaguchi, H.; Voiry, D.; Fujita, T.; Chen, M.; Chhowalla, M. *Nano Lett.* **2011**, *11*, 5111–5116. doi:10.1021/nl201874w

18. Ribeiro-Soares, J.; Almeida, R. M.; Barros, E. B.; Araujo, P. T.; Dresselhaus, M. S.; Cançado, L. G.; Jorio, A. *Phys. Rev. B* **2014**, *90*, 115438. doi:10.1103/physrevb.90.115438

19. Wilson, J. A.; Yoffe, A. D. *Adv. Phys.* **1969**, *18*, 193–335.
doi:10.1080/00018736900101307

20. Lin, M.-L.; Tan, Q.-H.; Wu, J.-B.; Chen, X.-S.; Wang, J.-H.; Pan, Y.-H.; Zhang, X.; Cong, X.; Zhang, J.; Ji, W.; Hu, P.-A.; Liu, K.-H.; Tan, P.-H. *ACS Nano* **2018**, *12*, 8770–8780. doi:10.1021/acsnano.8b05006

21. Cao, Y.; Fatemi, V.; Fang, S.; Watanabe, K.; Taniguchi, T.; Kaxiras, E.; Jarillo-Herrero, P. *Nature* **2018**, *556*, 43–50. doi:10.1038/nature26160

22. Yankowitz, M.; Chen, S.; Polshyn, H.; Zhang, Y.; Watanabe, K.; Taniguchi, T.; Graf, D.; Young, A. F.; Dean, C. R. *Science* **2019**, *363*, 1059–1064. doi:10.1126/science.aav1910

23. Autere, A.; Jussila, H.; Dai, Y.; Wang, Y.; Lipsanen, H.; Sun, Z. *Adv. Mater.* (Weinheim, Ger.) **2018**, *30*, 1705963.
doi:10.1002/adma.201705963

24. Hsu, W.-T.; Zhao, Z.-A.; Li, L.-J.; Chen, C.-H.; Chiu, M.-H.; Chang, P.-S.; Chou, Y.-C.; Chang, W.-H. *ACS Nano* **2014**, *8*, 2951–2958. doi:10.1021/nn500228r

25. Yu, H.; Liu, G.-B.; Tang, J.; Xu, X.; Yao, W. *Sci. Adv.* **2017**, *3*, e1701696. doi:10.1126/sciadv.1701696

26. Zhang, X.; Qiao, X.-F.; Shi, W.; Wu, J.-B.; Jiang, D.-S.; Tan, P.-H. *Chem. Soc. Rev.* **2015**, *44*, 2757–2785. doi:10.1039/c4cs00282b

27. Astié, V.; Wasem-Klein, F.; Makhlouf, H.; Paillet, M.; Huntzinger, J.-R.; Sauvajol, J.-L.; Zahab, A.-A.; Juillaguet, S.; Contreras, S.; Voiry, D.; Landois, P., in press.

28. Li, S.-L.; Miyazaki, H.; Song, H.; Kuramochi, H.; Nakaharai, S.; Tsukagoshi, K. *ACS Nano* **2012**, *6*, 7381–7388.
doi:10.1021/nn3025173

29. Lee, C.; Yan, H.; Brus, L. E.; Heinz, T. F.; Hone, J.; Ryu, S. *ACS Nano* **2010**, *4*, 2695–2700. doi:10.1021/nn1003937

30. Molina-Sanchez, A.; Wirtz, L. *Phys. Rev. B* **2011**, *84*, 155413.
doi:10.1103/physrevb.84.155413

31. Li, X.-L.; Qiao, X.-F.; Han, W.-P.; Zhang, X.; Tan, Q.-H.; Chen, T.; Tan, P.-H. *Nanotechnology* **2016**, *27*, 145704.
doi:10.1088/0957-4484/27/14/145704

32. Zhao, Y.; Luo, X.; Li, H.; Zhang, J.; Araujo, P. T.; Gan, C. K.; Wu, J.; Zhang, H.; Quek, S. Y.; Dresselhaus, M. S.; Xiong, Q. *Nano Lett.* **2013**, *13*, 1007–1015. doi:10.1021/nl304169w

33. Zhang, X.; Han, W. P.; Wu, J. B.; Milana, S.; Lu, Y.; Li, Q. Q.; Ferrari, A. C.; Tan, P. H. *Phys. Rev. B* **2013**, *87*, 115413.
doi:10.1103/physrevb.87.115413

34. Lee, J.-U.; Cheong, H. *J. Raman Spectrosc.* **2018**, *49*, 66–75.
doi:10.1002/jrs.5200

35. Dragcho, S. G.; Kim, J.; Hong, X.; Jin, C.; Shi, S.; Tongay, S.; Wu, J.; Wang, F. *Phys. Rev. B* **2017**, *95*, 165417.
doi:10.1103/physrevb.95.165417

36. Carvalho, B. R.; Malard, L. M.; Alves, J. M.; Fantini, C.; Pimenta, M. A. *Phys. Rev. Lett.* **2016**, *116*, 089904.
doi:10.1103/physrevlett.116.089904

37. Nemanich, R. J.; Tsai, C. C.; Connell, G. A. N. *Phys. Rev. Lett.* **1980**, *44*, 273–276. doi:10.1103/physrevlett.44.273

38. Lui, C. H.; Ye, Z.; Ji, C.; Chiu, K.-C.; Chou, C.-T.; Andersen, T. I.; Means-Shively, C.; Anderson, H.; Wu, J.-M.; Kidd, T.; Lee, Y.-H.; He, R. *Phys. Rev. B* **2015**, *91*, 165403.
doi:10.1103/physrevb.91.165403

39. Huang, S.; Liang, L.; Ling, X.; Puretzky, A. A.; Geohegan, D. B.; Sumpter, B. G.; Kong, J.; Meunier, V.; Dresselhaus, M. S. *Nano Lett.* **2016**, *16*, 1435–1444. doi:10.1021/acs.nanolett.5b05015

40. Liao, M.; Wei, Z.; Du, L.; Wang, Q.; Tang, J.; Yu, H.; Wu, F.; Zhao, J.; Xu, X.; Han, B.; Liu, K.; Gao, P.; Polcar, T.; Sun, Z.; Shi, D.; Yang, R.; Zhang, G. *Nat. Commun.* **2020**, *11*, 2153.
doi:10.1038/s41467-020-16056-4

41. Quan, J.; Linhart, L.; Lin, M.-L.; Lee, D.; Zhu, J.; Wang, C.-Y.; Hsu, W.-T.; Choi, J.; Embley, J.; Young, C.; Taniguchi, T.; Watanabe, K.; Shih, C.-K.; Lai, K.; MacDonald, A. H.; Tan, P.-H.; Libisch, F.; Li, X. *Nat. Mater.* **2021**, *20*, 1100–1105.
doi:10.1038/s41563-021-00960-1

42. Tiberj, A.; Rubio-Roy, M.; Paillet, M.; Huntzinger, J.-R.; Landois, P.; Mikolasek, M.; Contreras, S.; Sauvajol, J.-L.; Dujardin, E.; Zahab, A.-A. *Sci. Rep.* **2013**, *3*, 2355. doi:10.1038/srep02355

43. Yan, R.; Simpson, J. R.; Bertolazzi, S.; Brivio, J.; Watson, M.; Wu, X.; Kis, A.; Luo, T.; Hight Walker, A. R.; Xing, H. G. *ACS Nano* **2014**, *8*, 986–993. doi:10.1021/nn405826k

44. Melnikova-Kominkova, Z.; Jurkova, K.; Vales, V.; Drogowska-Horná, K.; Frank, O.; Kalbac, M. *Phys. Chem. Chem. Phys.* **2019**, *21*, 25700–25706. doi:10.1039/c9cp04993b

45. Lee, T.; Choi, J.-H.; Ahn, J.-H.; Yoon, Y.-G.; Rho, H. *Appl. Surf. Sci.* **2022**, *579*, 152208. doi:10.1016/j.apsusc.2021.152208
46. Liu, K.; Zhang, L.; Cao, T.; Jin, C.; Qiu, D.; Zhou, Q.; Zettl, A.; Yang, P.; Louie, S. G.; Wang, F. *Nat. Commun.* **2014**, *5*, 4966.
doi:10.1038/ncomms5966

47. Rice, C.; Young, R. J.; Zan, R.; Bangert, U.; Wolverson, D.; Georgiou, T.; Jalil, R.; Novoselov, K. S. *Phys. Rev. B* **2013**, *87*, 081307. doi:10.1103/physrevb.87.081307

48. Li, Z.; Lv, Y.; Ren, L.; Li, J.; Kong, L.; Zeng, Y.; Tao, Q.; Wu, R.; Ma, H.; Zhao, B.; Wang, D.; Dang, W.; Chen, K.; Liao, L.; Duan, X.; Duan, X.; Liu, Y. *Nat. Commun.* **2020**, *11*, 1151.
doi:10.1038/s41467-020-15023-3

49. Cortijo-Campos, S.; Prieto, C.; De Andrés, A. *Nanomaterials* **2022**, *12*, 1330. doi:10.3390/nano12081330

50. Liu, Q.; Li, L.; Li, Y.; Gao, Z.; Chen, Z.; Lu, J. *J. Phys. Chem. C* **2012**, *116*, 21556–21562. doi:10.1021/jp307124d

51. Debnath, R.; Maity, I.; Biswas, R.; Raghunathan, V.; Jain, M.; Ghosh, A. *Nanoscale* **2020**, *12*, 17272–17280. doi:10.1039/c9nr09897f

52. Mignuzzi, S.; Pollard, A. J.; Bonini, N.; Brennan, B.; Gilmore, I. S.; Pimenta, M. A.; Richards, D.; Roy, D. *Phys. Rev. B* **2015**, *91*, 195411. doi:10.1103/physrevb.91.195411

53. Zhou, X.; Jin, K.; Cong, X.; Tan, Q.; Li, J.; Liu, D.; Luo, J. *J. Colloid Interface Sci.* **2019**, *538*, 159–164.
doi:10.1016/j.jcis.2018.11.032

54. Maino, J. L.; Schouten, R.; Umina, P. *J. Appl. Ecol.* **2021**, *58*, 789–800.
doi:10.1111/1365-2664.13812

## License and Terms

This is an open access article licensed under the terms of the Beilstein-Institut Open Access License Agreement (https://www.beilstein-journals.org/bjnano/terms), which is identical to the Creative Commons Attribution 4.0 International License (https://creativecommons.org/licenses/by/4.0). The reuse of material under this license requires that the author(s), source and license are credited. Third-party material in this article could be subject to other licenses (typically indicated in the credit line), and in this case, users are required to obtain permission from the license holder to reuse the material.

The definitive version of this article is the electronic one which can be found at:
https://doi.org/10.3762/bjnano.15.26