Approach to determine measurement uncertainty in complex nanosystems with multiparametric dependencies and multivariate output quantities

In many cases, the determination of the measurement uncertainty of complex nanosystems provides unexpected challenges. This is in particular true for complex systems with many degrees of freedom, i.e. nanosystems with multiparametric dependencies and multivariate output quantities. The aim of this paper is to address specific questions arising during the uncertainty calculation of such systems. This includes the division of the measurement system into subsystems and the distinction between systematic and statistical influences. We demonstrate that, even if the physical systems under investigation are very different, the corresponding uncertainty calculation can always be realized in a similar manner. This is exemplarily shown in detail for two experiments, namely magnetic nanosensors and ultrafast electro-optical sampling of complex time-domain signals. For these examples the approach for uncertainty calculation following the guide to the expression of uncertainty in measurement (GUM) is explained, in which correlations between multivariate output quantities are captured. To illustate the versatility of the proposed approach, its application to other experiments, namely nanometrological instruments for terahertz microscopy, dimensional scanning probe microscopy, and measurement of concentration of molecules using surface enhanced Raman scattering, is shortly discussed in the appendix. We believe that the proposed approach provides a simple but comprehensive orientation for uncertainty calculation in the discussed measurement scenarios and can also be applied to similar or related situations.

In many cases, the determination of the measurement uncertainty of complex nanosystems provides unexpected challenges. This is in particular true for complex systems with many degrees of freedom, i.e. nanosystems with multiparametric dependencies and multivariate output quantities. The aim of this paper is to address specific questions arising during the uncertainty calculation of such systems. This includes the division of the measurement system into subsystems and the distinction between systematic and statistical influences. We demonstrate that, even if the physical systems under investigation are very different, the corresponding uncertainty calculation can always be realized in a similar manner. This is exemplarily shown in detail for two experiments, namely magnetic nanosensors and ultrafast electro-optical sampling of complex time-domain signals. For these examples the approach for uncertainty calculation following the guide to the expression of uncertainty in measurement (GUM) is explained, in which correlations between multivariate output quantities are captured. To illustate the versatility of the proposed approach, its application to other experiments, namely nanometrological instruments for terahertz microscopy, dimensional scanning probe microscopy, and measurement of concentration of molecules using surface enhanced Raman scattering, is shortly discussed in the appendix. We believe that the proposed approach provides a simple but comprehensive orientation for uncertainty calculation in the discussed measurement scenarios and can also be applied to similar or related situations.
Original content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
The Guide to the Expression of Uncertainty in Measurement (GUM) was published in 1995 to provide a well defined and mathematically sound procedure for the determination of measurement uncertainty as an essential part of every measurement result in general [1]. It has become the standard procedure for all measurements performed in metrology institutes and accredited calibration laboratories and is increasingly in use in research institutes, universities and industry. It serves as well as the basis for accurate process control and thus quality management in industrial production. In the meantime it has been extended by two supplements to clarify and extend its methodology for situations where Monte Carlo simulations are applied for the propagation of probability distributions (Supplement 1) [2] or a larger number of output quantities is dealt with (Supplement 2) [3].
Although the GUM and its supplements cover many experimental scenarios, uncertainty propagation for complex nanosystems is still challenging. This is mainly because the measurement system, consisting of the instruments and the investigated nanosystem, provide a huge number of input and output parameters and are very difficult to model [4]. Moreover, correlations in such multivariate problems might significantly influence uncertainty propagation and, thus, need to be considered during data analysis.
Yet, complex multiparametric measurement scenarios where correlations play an important role can often be simplified such that uncertainty propagation can be performed with an acceptable effort. This will be demonstrated in this paper. This work is the result of discussions in the Research Training Group GrK1952/1 'Metrology for Complex Nanosystems' of the Deutsche Forschungsgemeinschaft (DFG), where all participating projects deal with the determination of the measurement uncertainty of complex nanosystems.
The remainder of this paper is structured as follows. In section 2 we outline the methodical approach which includes the discussion of the modified Ishikawa diagrams and the propagation of uncertainties following the GUM. After that, we describe a detailed application of the aforementioned method for magnetic nanosensors and ultrafast electro-optical sampling of complex time-domain signals in section 3. This includes a full analysis of the uncertainty propagation. Finally, conclusions are given in section 4. To explain the presented approach in more detail, we discuss three additional applications (nanometrological instruments for THz microscopy, dimensional scanning probe microscopy and measurement of concentration of molecules using surface enhanced Raman scattering) in the appendix, where we schematically describe important parameters for uncertainty propagation.

Methodic approach for uncertainty calculation of complex systems
Conventional methods to determine measurement uncertainty may only be applied to a certain set of problems. If a measurement system reaches a certain degree of complexity, the methods need to be modified or a completely different approach must be chosen. In order to be able to judge the complexity of a system a number of measures and criteria will be introduced here.
As a very basic definition, a complex system is something that is neither easy to model, nor easy to manufacture. This places it somewhere between completely ordered (easy to model, possibly difficult to manufacture) and completely disordered (easy to manufacture, difficult to model adequately) [5]. General criteria that are also relevant in the determination of measurement uncertainty are described in more detail in the literature [6]: • A property requiring a novel approach is the nonlinearity of a system. Nonlinearity means that the superposition principle for input and output quantities on which the classical approach of the GUM relies is no longer valid. If a nonlinear system can not be sufficiently linearized in its operating point, a Monte Carlo approach must be chosen [2]. • Another criterion is the consideration of feedback mechanisms, which means that output quantities not only rely on input quantities, but also on output quantities of previous states of the system. In this case the dependence of the output parameters on the input parameters of a system cannot be described in an analytical way. • A related property of complex systems which can complicate uncertainty calculations, is called emergence and essentially describes complex properties that a system exhibits caused by a high number of physically describable processes. An example of this is magnetic hysteresis resulting from the behaviour of magnetic domains. While a model, which gives a sufficient accuracy for the desription of the magnetization of a macroscale (ferro-)magnet exists, the microscale structure of magnetic domains is never the same during different magnetization and demagnetization cycles of a given material.

Modified Ishikawa diagrams
Ishikawa diagrams were created in 1968 by Kaoru Ishikawa, a Japanese organizational theorist, as cause-and-effect diagrams which describe the causes of a specific event [7]. They are commonly applied in product design and quality defect prevention to identify potential factors causing an overall effect.
The standard ISO 9004-4 gives the following description: 'A cause-and-effect diagram is used to analyze cause-and-effect relationships; to communicate cause-and-effect relationships; and to facilitate problem solving from symptom to cause to solution' [8].
Ishikawa diagrams start with a horizontal arrow, the main bone. Major categories of possible causes are represented using new arrows pointing to the main bone to identify the influences of the causes on the event. Factors to be considered in this process include data and information systems, environ ment, equipment, materials, measurements, methods, and observer [8].
From a metrology point of view, Ishikawa diagrams are a compact way to identify, illustrate and structure the effects of input quantities on an output quantity, which usually is the measurand of interest, and then facilitate the creation of the uncertainty budget for the measurement result.. Also the expression of system equations in mathematical form is facilitated. However, for a complex system, an Ishikawa diagram may have a large number of causes. Many contributing input parameters exist, for which the uncertainty propagation is typically very complicated. In order to simplify such systems and to gain a better understanding of such systems, a breakdown of the system into independent subsystems is desirable.
With the separation of the block diagram (figure 1) into systematic influences (subsystems above the main bone) and statistical influences (subsystems below the main bone), uncertainty propagation can be performed in a straightforward way. Error sources can be distinguished by their observed or expected behaviour. While statistic errors behave in a nondeterministic way and can not be predicted, systematic errors are expected to be the same for a specific set of input conditions. Statistical influences can be identified through repeated measurements, systematic errors through comparison measurements with different setups or devices. The systematic influences are captured by the system equation (see section 2.2) and the statistical influences are captured by assigning a best estimate and an uncertainty (or probability density function) to the input quantities. The uncertainty can either be obtained from repeated measurements (Type A uncertainty) or from scientific judgment based on other information, such as manufacturer's specifications or calibration certificates (Type B uncertainty) [1].

System equation
To propagate uncertainties from input variables to an output quantity, an equation for the system is needed. This system equation is also referred to as measurement model and describes the influence of all uncertainties of the input variables on the output quantity in an analytical way. For most systems this equation can be derived from the Ishikawa diagram. Yet, for complex systems this approach cannot be made intuitively because of the many parameters which have to be taken into account.
In the approach described in section 2.1 it is necessary to modify the Ishikawa diagram in a way that its subsystems are independent of each other. Under this condition it is possible to calculate the uncertainty of every subsystem h i independently and combine these uncertainties to the uncertainty of the main output quantity Y being related to the input variables X i through the system equation Hereby, f describes the system equation as a function of the input quantities, while f is a function of the uncorrelated subsystems (see figure 1). The division of the whole system into independent subsystems makes the calculation less complex and less time consuming. Yet, it still captures the main influencing factors and preserves the calculation of correlations in the output quantity since the best estimate of Y and its covariance matrix will be obtained. The calculation can be executed in two ways. Both of these techniques will be briefly described in the following.

Propagation of measurement uncertainty using covariance matrices
The algebraic method using covariance matrices poses some demands on the model [1]. The first condition is that the model is either linear or can be linearized in the working point using a Taylor series. Usually, the latter condition requires small uncertainties for the input parameters. The other requirement focuses on the probability density function (PDF). This function is restricted to symmetric and common PDFs [1].
The best estimate of Y is obtained from (1) using the best estimates of the input variables. The covariance matrix of Y is given by with the diagonal and the off-diagonal elements denoting the variances u(Y i , Y i ) and covariances u(Y i , Y j ), respectively, of the variables of Y. The covariance matrix can be calculated with the following generalized expression [3] where C X is the sensitivity matrix with the dimension (N × M) and is given by The covariance matrix of the input variables U X will be obtained from repeated measurements of the input variables.
In this way the statistical influences on the measurement process will be captured. For the special case of non-correlated subsystems, U X can be separated into a block-diagonal form, where the diagonal matrices correspond to the covariance matrices of the subsystems: An example for uncertainty propagation using this technique will be given in section 3.1. Yet, to better visualize correlations we will focus on the correlation matrix R Y rather than the covariance matrix U Y . The indices of the correlation matrix are obtained from with σ denoting the standard deviation of Y. The values of the correlation matrix always lie in the interval [−1,1]. For multivariate output quantities the specification of coverage intervals is not as easy as for a single scalar output quantity [3]. This is because a region R Y in the N-dimensional space of the output quantity Y is required such that the probability that Y lies in R Y is equal to the specified coverage interval. Once the PDF of Y is known, the coverage probability for a specified coverage interval can be obtained. This is straightforward for a Gaussian PDF, but numerical methods such as Monte Carlo calculations are required for other PDFs [3].

Propagation of measurement uncertainty using Monte Carlo simulations
In general, the Monte Carlo method describes a numerical procedure to solve equations by utilizing large numbers of computer generated random samples. In the field of uncertainty analysis, the Monte Carlo method is applied to calculate the probability density function of model systems (g Y ) from the probability density functions (g X1 , . . . , g XM ) of the associated input parameters, see figure 2. The advantage of the Monte Carlo method is that it is applicable for arbitrary input probability distributions even in cases of nonlinear system equations.
The procedure to employ the Monte Carlo method is described in GUM Supplement 1 [2] and extended in GUM Supplement 2 [3] for an arbitrary number of output parameters including a full correlation treatment [9,10] and can be summarized in the following five steps: The disadvantage of this method is that it requires extensive number of evaluations of the system equation (e.g. GUM Supplement 1 advised 200 000 Monte Carlo trials for 95% confidence). Especially for cases which include numerical expensive algorithms or for complex systems with many correlated input quantities, Monte Carlo simulations will not easily converge after a practicable period of time. One way to overcome this is by utilizing more adaptive methods [11] to determine the required number of Monte Carlo trials and effectively reduce the number of Monte Carlo runs.
The expanded uncertainties are typically specified by the corresponding coverage intervals (i.e. 95% coverage interval) derived directly from the PDF of the output quantity.

Applications
In this section, we demonstrate two detailed applications for the methods which were presented in the previous section. The starting point for both applications is the Ishikawa diagram which describes the influences of all input quantities on the multivariate output quantities. We explicitly discuss the simplification of the Ishikawa diagrams to block diagrams which have a considerably reduced complexity and are separated Figure 2. Scheme for the propagation of the probability density function. Here, g X1 , . . . , g XM and g Y1 , . . . , g YN correspond to the probability density functions of the input quantities and the output quantities, respectively. into systematic and statistical influences. From this, we calculate the best estimate, its uncertainty and the associated correlation matrix of multivariate output quantities for the two methods. Hereby, the required system equations to account for systematic influences are derived from the corresponding block diagrams. In the first subsection, we demonstrate the method for a magnetic nanosensor as an example of a sensor in frequency-domain. We linearize the derived system equation in order to directly utilize equation (4) for the propagation of the covariance matrix. Second, we illustrate results for ultrafast electro-optical sampling of complex time-domain signals using Monte Carlo simulations. In this case, no linearization of the system equation is required but the multivariate PDF of the input quantities, which might include correlations, have to be known. In addition, three more applications for different sensors are discussed in the appendix.
In the following, lower-and upper-case variables denote time-and frequency domain signals, respectively, with the time and frequency dependence being taken as implicit. This differs from the notation used in section 2, where only uppercase variables have been used.

Magnetic nanosensor
The magnetic noise versus frequency of a magnetic nanosensor is analyzed as the output quantity, which is also referred to as magnetic flux density noise. Magnetic nanosensor is a term for miniaturized magnetoresistive (MR) sensors. MR sensors, like anisotropic magnetoresistance (AMR), giant magnetoresist ance (GMR) and tunneling magnetoresistance (TMR) are types of magnetic field sensors that are widely used and commercially available. This is because of their low price, small dimension and-compared to hall sensors-the high signal-to-noise ratio [12,13]. In spite of that, the noise characteristics of these sensors and also the dependence of the noise characteristics on the manufacturing process are not completely known. Especially the 1/f-noise characteristics are rather unexplored, but very important, because the magnetic noise S 1/2 B determines the detection limit [12]. The measurement unit is T/ √ Hz. It is a combination of a sensitivity measurement (V/T) and a voltage noise measurement (V/ √ Hz). A high sensitivity and a low intrinsic voltage noise are desirable and lead to a high signal-to-noise ratio. In particular, the 1/f noise of GMR and TMR sensors are severe limitations, which influence the application significantly.
A very sensitive characterization system for AMR, GMR, and TMR sensors was developed in [14]. It facilitates the opportunity of noise measurements of different temperatures and magnetic bias fields. Thus it is possible to increase the knowledge about 1/f noise sources in MR sensors. But noise measurements themselves are a very complicated and critical measurement task. It is not possible to measure the intrinsic noise alone. The noise measurement always contains a combination of different noise sources: the noise of the sensor itself, the noise of the measurement electronics and environmental influences. It should be ensured that the intrinsic sensor noise is much higher than the other influences. Therefore, much experience about noise and intrinsic noise sources in electronics is necessary. Additionally, the sensor noise depends on the working point and the condition of the sensor. Parameters like temperature, thin film magnetization and bias field have to be taken into account and to be controlled very carefully. Typically, the sensor behaviour is nonlinear. In the working point, the sensor behaviour can be linearized and many nonlinear dependencies are negligible.
Also, the 1/f noise measurement result is difficult to evaluate. It is not generally possible to predict the influences of the manufacturing on the result with magnetic and electric simulations. In addition, there is no 1/f noise standard available, which could be used at low frequencies to quantify the results of the measurement system. Therefore, it is a very important task to elaborate the measurement uncertainty of magnetic noise measurements with the aid of the GUM. The measurements are influenced by many other factors like temperature, offset, input noise sources etc [14]. Influencing quantities of the corresponding measurement are collected in an Ishikawa diagram, see figure 3(a). In this figure it is shown that the magnetic noise is a combination of a sensitivity measurement and a noise measurement. These measurements are separated in subsystems and are evaluated independently. This is possible, because the noise measurement and the sensitivity measurement are analyzed in the same working point (bias field, temperature) and the sensitivity behavior is assumed to be linear in this working point. The sensitivity measurement can be separated into the magnetic field generation and the sensor output measurement. These systems are also independent.
The noise measurement can be separated into four subsystems. The external spectrum analyzer (Agilent 35670A) is assumed to be negligible because of the high input impedance, the intrinsic calibration system and the dynamic range adjustment. The influence of the environment is suppressed significantly with the aid of the temperature and shielding system. It can also be modeled as an independent system. A problem is the interaction between the sensor and the first amplifier stage. The intrinsic noise sources of the amplifier (voltage and current noise) superimpose the intrinsic noise of the sensor. In spite of that, the systems are separated into different and independent subsystems. This is possible, because there are good models to distinguish the influences very accurately [15].
Our formulated method presents not only the separation of the whole system into different subsystems, but also the Ishikawa diagram is simplified to a smaller diagram (section 2.1) which is used to derive the system equation. Therefore, a first approach can be the determination of the largest systematic and statistical influences.
These quantities are shown in the block diagram in figure 3(b). At first view, it seems curious that noise is arranged at systematic influences in the upper area of the diagram. Here, the systematic frequency characteristics of noise is meant, which is contained in the measurement result. It can be measured independently of the sensor sensitivity because of the previously explained measurement condition. Additionally, the intrinsic noise sources of the amplifier are systematic influences with a strong contribution to the output quantity. These different input quantities are used to deduce the output equation. The statistical influences are collected in the lower side of the diagram. These influences are difficult to capture and to model. The influences of these quantities were suppressed with the aid of the measurement system. The remaining uncertainties are captured by repeated measurements.
The output quantity of interest is the magnetic noise versus frequency, which is also referred to as magnetic flux density noise S 1/2 B and is given by The system equation is deduced from the block diagramm ( figure 3(b)). Here, S meas,V is the voltage noise of the nanosensor and the scalar S Sens denotes a frequency independent sensitivity factor. The amplifier used to amplify the nanosensor's signal has two noise contributions which need to be accounted for S Amp,V is the input voltage noise and S Amp,RI is the noise resulting from the input current. The frequency dependent values are taken from the datasheet [15] (Type B).

Uncertainty propagation for the magnetic nanosensor
We now present a way to calculate the best estimate and the corresponding uncertainty for magnetic nanosensors utilizing the covariance matrix formalism introduced in section 2.3. The measurements were performed on an AMR sensor of Sensitec GmbH [16]. The best estimate of S meas,V and its covariance matrix were obtained from repeated noise measurements in the frequency range from 10 Hz to 3.2 kHz. Thirty noise traces (each averaged 100 times) and 30 sensitivity measurements were recorded and analyzed. Therefore, the corresponding uncertainty contributions are Type A contributions. This data was merged in the covariance matrix of the subsystems using equation (6).
The semitransparent color in figure 4(a) denotes the best estimate ± the expanded uncertainty [1] with a coverage factor k = 2 such that a 95% confidence level is obtained. The white noise characteristic is observable in the higher frequency region. In the lower frequency band, the 1/f noise dominates and a frequency dependent characteristic of the sensor noise is obtained. It is known that the uncertainty of voltage noise measurements increases for lower frequencies which is caused by the time scale. For low frequencies, only a few values contribute to the result. Usually, a further measurement with a smaller frequency limit would be recorded and a merged spectrum would be shown. Then frequencies smaller than 100 Hz would be replaced with a more representative measurement. We decided to renounce on this since this is not required for the general discussion of noise measurements.
The covariance matrix was calculated from (6) and the covariance matrix of the output quantity was obtained with equation (3). This matrix was normalized (correlation coefficient) by a division with the product of the standard deviation. Here, the correlation coefficient matrix is a good concept to evaluate the noise measurement ( figure 4(b)). The white noise which is located in the high frequency range is not correlated. The 1/f noise is dominating in the small frequency range. The correlation of 1/f noise is still an important research subject and not completely understood. But the mathematical model and the calculated cross correlation show a nonlinear and decayed behaviour [17]. So, a linear correlation would not be expected. In contrast, disturbing signals like drift signals and electromagnetic signals may be correlated. Altogether, the correlation coefficient matrix allows further evaluation of the result.
The correlations between the different frequency points are visualized by the correlation matrix in figure 4(b). While the frequency scale of the matrix is logarithmic, the color scale itself is linear. Although the diagonal elements of the correlation matrix are always one, this relationship is not apparent from figure 4(b) due to the finite pixel width of the plot and the logarithmic frequency scale. Significant correlations are visible in the off-diagonal elements of the matrix with absolute values of up to 0.4 and a small systematic pattern is visible at low frequencies (where 1/f noise dominates). This is not caused by the correlation of the intrinsic noise of the sensor. It is a measurement artifact which is caused by the changing DC-level of the output signal in connection with the limited number of measurement points in the low frequency range. It is another indication for the complexity of noise measurement. Further measurements on larger time scales are necessary to replace this area as explained before. At higher frequencies (where white noise dominates) the correlation between the different frequency points is stochastically distributed.

Ultrafast electro-optical sampling of complex signals in time domain
The laser-based sampling of complex time-domain voltage signals, representing the output quantity under study, is realized by employing a femtosecond laser and a so-called pump-probe technique. With this technique the voltage signals propagating on transmission lines such as coplanar waveguides and comprising several thousand data points can be detected with a bandwidth exceeding 1 THz [18,19]. Uncertainty propagation for these signals requires elaborated signal processing algorithms [20]. In particular, the accurate consideration of correlations between the several thousand data points is a challenging task [21,22]. This is obvious since a correlation matrix of 1000 input variables has 1 000 000 data points and several million measurements are required for a full-rank correlation matrix.
While the generation and detection of ultrashort voltage pulses can be accomplished without the involvement of nanostructures, recent work focuses on the usage of semiconductor nanostructures for the quantitative detection of electric fields. This establishes the link between electro-optic sampling (EOS) and nanosystems.
The Ishikawa diagram of the electro-optic sampling technique is shown in figure 5(a). The output quantity is influenced by systematic and statistic effects with three main causes: the femtosecond laser, optical components, and electronics. Starting from this Ishikawa diagram we separate influencing factors contributing to the time trace of an ultra-short voltage pulse into systematic and statistical influences. Both groups can be divided into uncorrelated subsystems. This simplified Ishikawa diagram is presented in the block diagram figure 5(b).
The stochastic uncertainty contributions may be quantified by multiple measurements while keeping the experimental configuration constant. The main systematic influences are the measured voltage pulses. This subsystem includes several measurements as noted below. Independent from that, we can identify a systematic error which arises from a varying time base and a varying intensity of the measurement signal for different sampling times. This influence forms the subsystem EOS-reference signal. The last subsystem (EOS-transfer function) is an error contribution which results from light-substrate interaction and a finite pulse width of the probe beam acting as a low-pass filter and, thus, leading to a broadening of the ultrashort voltage pulse. This process is described by a convolution with an electro-optic transfer function.
These systematic influences can be corrected, if their influence on the output quantity is known (although the uncertainty of the correction might have an influence on the output quantity). The associated system equation to remove the systematic influences and to calculate the time trace of an ultrashort voltage pulse v as derived from the block diagram in figure 5 A detailed derivation of the system equation and an explanation of the input quantities is given in [19]. Here, we only briefly discuss the input quantities and the relation to the subsystems in the block diagram. First, the measured voltage pulse v meas is normalized by the EOS-reference signal v ref in time domain.
After that the normalized pulse is deconvolved with the EOS-transfer function H EOS . This is done in frequencydomain by utilizing the Fourier transformation (F ). Finally, we have to account for multiple reflections of the measured voltage pulse. The calculation of the reflection coefficient Γ requires the measurement of several voltage pulses and is presented in detail in [19]. Again, the deconvolution of the normalized voltage pulse with the reflection coefficient is executed in frequency-domain.

Uncertainty propagation for ultrafast electro-optical sampling
To obtain the time trace of an ultra-short voltage pulse propagation on a coplanar waveguide the system equation (9) is employed. The corresponding uncertainty propagation is performed by the Monte Carlo simulation described in the section 2.4. The reason for using the Monte Carlo algorithm is that the system equation includes deconvolution operations utilizing the Fourier transformation. For this case a direct uncertainty propagation with the covariance matrix is more complex and demanding than the Monte Carlo approach.
The mean value and the standard deviation of the measured voltage pulse v meas is calculated from multiple measurements under constant experimental conditions. From this, the PDF is estimated using a t-distribution. The reference signal v ref is measured by applying a constant external bias to account for intensity changes of the optical probe beam for different sampling times. Since the temporal shape of v ref can be described using few parameters (see [19]) the complete covariance matrix can be estimated. Hence, the Monte Carlo realizations for v ref are drawn from a multivariate normal distribution.
To account for the electro-optical transfer function H EOS , we divide the normalized voltage pulse by the electro-optic transfer function in frequency-domain which acts as a deconvolution operation. Correlations between H EOS and the measured voltage pulse are neglected so that we can treat the uncertainty of the transfer function in an independent subsystem. Furthermore, we have to consider the multiple reflections of the measured voltage pulse with the help of the reflection coefficient Γ. The equation to calculate the reflection coefficient Γ (see [19]) uses voltage pulse measured at different positions on the coplanar waveguide. Therefore, the t-distributed Monte Carlo realizations from v meas are used to calculate the PDF of the reflection coefficient.
Uncertainty contributions arising from instruments like lock-in amplifiers and signal generators are much smaller than the uncertainties described above and are therefore neglected. Thus, all uncertainty contributions are obtained from repeated measurements and only Type A uncertainties appear in the uncertainty analysis in this subsection.
As a result, we receive a multivariate PDF for the voltage pulse vector v. The mean value and correlation matrix of these output quantities can be directly calculated from this multivariate PDF, see figure 6. The best estimate in figure 6(a) Figure 6. (a) Best estimate of the calculated ultrashort voltage pulse v according to equation (9). Systematic influences which alter the shape of the pulse have been removed. The semitransparent color denotes the 95% confidence interval. This interval is so small that it only becomes visible at the peak of the voltage pulse. (b) Corresponding correlation matrix. displays an ultrashort voltage peak at 25 ps with a pulse width of approximately 2 ps. The semitransparent color in figure 6(a) denotes the 95% confidence interval as directly obtained from the Monte Carlo calculations.
The corresponding correlation matrix is shown in figure 6(b). Again, we focus on the correlation matrix rather than on the covariance matrix to better visualize correlations. Significant correlation values can be observed along the axes at 25 ps, emphasizing the influence of the main peak on voltage values at other time instances. This behavior can also be observed for smaller voltage pulses for example at 105 ps resulting in a certain pattern of the correlation matrix.

Conclusions
In order to find a best estimate and the corresponding uncertainty of a multivariate output quantity, it is often useful to reduce the level of complexity of the system, while simultaneously considering correlations in input and output quantities.
Here, we presented a new approach for the conversion of rather complex Ishikawa diagrams into simple block diagrams by dividing the system into subsystems and distinguishing between systematic and statistical influences. The block diagram allows for the uncertainty calculation according to GUM and captures correlations in the output quantity.
It should be emphasized that correlations in multivariate output quantities are always important if the output quantity serves as an input quantity for other calculations, i.e. in case of uncertainty propagation. Another equally important issue is the calculation of single parameters from a multivariate output quantity. In such a case, correlations might have a significant influence on the resulting single parameter. An example for such a case is the temporal width of the voltage pulses discussed in section 3.2. Correlations in the covariance matrix of the voltage pulse might significantly change the uncertainty of the temporal width.
Despite the specific examples given in this work, the presented guidelines will also be applicable to other laboratory setups and experiments, three of such additional examples are shortly discussed in the appendix. Thus, we are confident that our work will be helpful for other groups working on uncertainty evaluation of complex nanosystems for which correlations in multivariate input and output quantities occur.
The THz microscope [23,24] uses a Josephson cantilever [25] together with a positioning system to measure the spatial power distribution of high-frequency electromagnetic radiation up to the THz regime. The sensor is called a Josephson cantilever due to the fact that the Josephson effect of the hightemperature superconductor YBa 2 Cu 3 O 7 is used to generate unique IV-characteristics [26]. These IV characteristics are the output quantity for the determination of the measurement uncertainty. When high-frequency radiation couples into the cantilever, the IV characteristics change such that the power and the frequency of the radiation can be calculated [27].
An IV characteristic is recorded at each spatial position and three dimensional visualizations of the power distribution of the radiation can be generated [28]. Samples operated at high frequencies can be examined as well as passive elements. For this, a FIR laser can be used as THz source [29].
An overview of the influencing factors on the measurement uncertainty of the IV characteristics is given in the Ishikawa diagram in figure A1(a). The superconducting Josephson cantilever is not thoroughly describable by a single model as the effect of superconductivity for high temperature superconductors is not fully understood and many factors can have an influence on the IV characteristics. The Josephson canti lever is for example very sensitive to magnetic fields and temperature changes. The properties of the YBa 2 Cu 3 O 7 thin layer vary due to the production process and due to aging. To enable traceable measurements with the setup, a deeper understanding and a quantification of the influencing factors is necessary.
The effects on the IV characteristic of many influencing factors can be described in the RCSJ model [30,31]. A simulation of these effects is possible by numerically solving the according differential equations. The effects can be separated into subsystems to simplify the calculations. The most relevant subsystems are shown in figure A1(b). These subsystems are separated into systematic and statistical influences as proposed in chapter 2.1.

A.2. Dimensional STM-nanometrology based on atomic lattice constants
A metrological STM, which provides traceable measurements of sub-nanometer lengths and step heights is currently under construction. The projected field of application for this instrument is the undergoing development of new length standards in the nanometer regime, which are based on principles of self-assembly and invariants of nature. One of such invariant of nature is the 7 × 7 reconstruction of the Si(1 1 1) unit cell. With a lattice constant in the range of 5 nm [32] it is ideally suited as the base for calibration standards for scanning probe microscopes that cannot reach atomic resolution. In order to provide measurements of the aforementioned lattice constant as well as a method to control the manufacturing quality of newly developed nanostandards, an ultra-high vacuum scanning tunneling microscope (STM) was upgraded with interferometers, which cover all six degrees of freedom.
The STM is an instrument that can easily achieve atomic resolution [33], while the interferometers provide direct traceability of the length measurements to the international system of units (SI) [34]. The interferometers are based on a heterodyne interferometer design already existing and in use at the PTB, which has been modified with regard to operation under UHV conditions as well as the limited space available in the microscope assembly.
Sources of measurement uncertainty in the instrument are shown in the Ishikawa diagram, see figure A2(a). They consist of errors caused by interferometric length measurement together with those caused by the electronic devices and comp onents used to analyze the output signals of the interferometers [35].
Other errors such as the cosine and Abbe errors are inherent to length measurement and can be influenced and magnified by instabilities and tolerances of the mechanical assembly. A first separation of major error sources into systematic and statistic influences is illustrated in the block diagram ( figure A2(b)) and can be used as a basis for an analysis of uncertainty.
Certain subsystems, most prominently the interferometers [36], exhibit nonlinear behavior and can thus not be analyzed by applying classical methods, as described in section 2. Separate analysis of the relevant subsystems by utilizing a Monte Carlo method can be a solution. Furthermore, setting up a complete system equation for a scanning probe microscope is not feasible, making alternative methods for uncertainty calculation necessary [37].

A.3. Measurement of concentration of molecules by isotope dilution surface-enhanced Raman scattering
Generally speaking, in an isotope dilution surface-enhanced Raman scattering technique, intensity ratio of Raman peaks between target analyte and its isotopologue (as an internal standard) will be measured. With that, the unknown concentration of the target analyte can be determined from the intensity ratio and the known concentration of the reference isotopologue.
Raman spectroscopy is an optical technique that is based on inelastic scattering of monochromatic light. The so-called Raman shift is used in chemical analysis to identify materials via their molecular vibrations ('fingerprint' of a sample). When a sample is mixed with a isotopologue which is identical to the target molecule but isotopically enriched in at least one isotope. Due to the isotopic enrichment a slight increase in molar mass of the spike molecule leads a spectral shift of vibrational modes towards lower wavenumbers. Therefore, two techniques, Raman spectroscopy (RS) and isotope dilution, are combined allowing a sensitive, traceable and accurate concentration measurement that could compete with the more expensive isotope dilution mass spectrometry (ID-MS) [4,40]. In ID-RS case, both sample and internal standard are nearly equally affected by the complex contributions to the measurement uncertainty so that these contributions are being nearly eliminated in analogy to ID-MS.  However, ordinary Raman scattering is a rather weak process compared to Rayleigh scattering and light intensity measurements are elaborate [38,39]. In contrast, in case of surface-enhanced Raman spectroscopy (SERS), metal (e.g. Au or Ag) nanostructures on a substrate provide locally enhanced electromagnetic fields, so-called hot-spots, being capable of amplifying the Raman response by several orders of magnitude [41]. After adding molecules onto such a SERS substrate, each molecule will gain an enhancement factor k and the intensity contribution from each molecule would be kΩ, where Ω is the ordinary Raman cross section of one molecule.
Tiny variations of parameters of one single nano unit in principle will nonlinearly fluctuate the plasmonic property of the whole system and finally can change the SERS response [42]. In addition, the sample preparation to isolate the target molecule from the surrounding matrix for SERS measurements, e.g. in clinical or environmental applications, has a tremendous influence on Raman and SERS scattering processes, respectively [43]. Thus, it is neither easy nowadays to exactly model such a many-component system, nor easy to manufacture ideal metallic structures at the sub-nanoscale in a wellcontrolled manner.
Therefore, the breakdown of the system as shown in figure A3(a) would be considered as a candidate approach to determine the efficiency and uncertainty related to SERS. With that we confine the most complexities into the 'SERS substrate' subsystem and deal with them for instance using the Monte Carlo method, clarifying the uncertainty analysis within the remaining subsystems [44]. In figure A3(b) the simplified block diagram splits the many-component effects from the Ishikawa diagram into measurement uncertainties related statistic effects and to those that are related to corrections for systematic effects.
The isotope dilution SERS (ID-SERS) can be used as an extended SERS technique for traceable concentration measurements of a target molecule normally in its natural isotopic composition and is based on the principle that an isotopically enriched analyte (spike) is used as an internal standard [40]. Hence, considering ID-SERS measurements, the intensity ratio of the isotopologues is used in the Raman spectra for quantifying the concentration of the target molecule in combination with multivariate data analysis [45]. The major contributors to the uncertainty result from the knowledge of the isotopically labelled spike material (gravimetrically determined amount of substance, purity, etc). These spike material related effects, in turn, can be further reduced by applying double isotope dilution procedures as described in [4] and lead to similar measurement uncertainties as achieved by ID-MS as shown in [40].