Asymmetrical Transport Distribution Function: Skewness as a Key to Enhance Thermoelectric Performance

How to achieve high thermoelectric ﬁ gure of merit is still a scienti ﬁ c challenge. By solving the Boltzmann transport equation, thermoelectric properties can be written as integrals of a single function, the transport distribution function (TDF). In this work, the shape e ﬀ ects of transport distribution function in various typical functional forms on thermoelectric properties of materials are systematically investigated. It is found that the asymmetry of TDF, characterized by skewness, can be used to describe universally the trend of thermoelectric properties. By de ﬁ ning symmetric and asymmetric TDF functions, a novel skewness is then constructed for thermoelectric applications. It is demonstrated, by comparison with ab initio calculations and experiments, that the proposed thermoelectric skewness not only perfectly captures the main feature of conventional skewness but also is able to predict the thermoelectric power accurately. This comparison con ﬁ rms the unique feature of our proposed thermoelectric skewness, as well as its special role of connection between the statistics of TDF and thermoelectric properties of materials. It is also found that the thermoelectric performance can be enhanced by increasing the asymmetry of TDF. Finally, it is also interesting to ﬁ nd that the thermoelectric transport properties based on typical quantum statistics (Fermi-Dirac distributions) can be well described by typical shape parameter (skewness) for classical statistics.


Introduction
Recently, the thermoelectric energy conversion technology attracts more and more attentions due to its ability to convert heat into electricity or vice versa, with extreme reliability and without producing greenhouse gas emissions [1][2][3][4][5][6]. Therefore, thermoelectric devices can be used for power generation or refrigeration devices [1][2][3][4][5][6][7][8]. The energy conversion efficiency of thermoelectric devices is evaluated in terms of a dimensionless thermoelectric figure of merit, ZT, which depends on the transport coefficients of the constituent materials, and can be defined by [6].
where σ is the electrical conductivity, S is the Seebeck coefficient (thermoelectric power), T is the mean operating temperature, κ e and κ l are electronic and lattice thermal conductivities, respectively. The term σS 2 , the so-called power factor, characterizes the pure electronic contribution to ZT. A general strategy for improving thermoelectric efficiency is to increase the electrical conductivity and Seebeck coefficient (equivalently, to increase the power factor) but to decrease the thermal conductivity [5]. However, in real materials, the improvement of thermoelectric efficiency is not straightforward but rather difficult, due to the strong coupling among three associated components, σ, S, and κ. For example, it is well known that, compared with insulator, metals have very good electrical conductivity and one may expect from Equation (1) that metals can be good thermoelectrics. In fact, most metals are very poor thermoelectric materials due to their very small Seebeck coefficient and large thermal conductivity. For insulators with large band gaps, they have very large Seebeck coefficients. But their ultra small electrical conductivity significantly reduces the power factor and thus results in a very small ZT [5]. The optimal thermoelectric materials are therefore located in the region near the crossover between the semiconductor and metal, and the mechanism for enhancement of thermoelectric performance is thus desirable to be further studied. Different approaches have been reported to increase the thermoelectric figure of merit. To reduce the lattice thermal conductivity, increased phonon scattering by defects and nanocrystals [9][10][11], by layer structured PbSnS 2 or doping [12] in PbTe-based thermoelectric materials [13], and by the hierarchical phonon scattering centers [14,15] and interfaces [16][17][18], as well as the nanowires [19], has been reported. The electronic conductivity, Seebeck coefficient, and power factors have been studied by investigating the effects of pressure [20,21], strains [22][23][24], defects [25], doping [12,26], polytypes [27], superlattices [28,29], and reduced dimension or nanostructures [30,31]. The effects of symmetry of crystal structure on thermoelectric properties [32,33] and the effects of Aubry-Andre-Harper modulation on nanoscale thermoelectrics [34] have been investigated. Advancements have been made recently for the Bi 2 Te 3 nanofilms [28,35], other chalcogenide-based materials [26,[36][37][38], CaMg 2 Bi 2 -based compounds [39] or Mg 3 Sb 2- x Bi x alloys [40,41], half-Heusler materials [42], semimetals [43], and organic materials [44][45][46][47]. An interesting increasing trend can be observed by utilizing the machine learning method to assist research in energy materials [48][49][50] including thermoelectric materials [42,51].
All the aforementioned approaches are aimed at optimizing some particular components of thermoelectric properties. Unfortunately, these components are coupled with each other and are not easy to be optimized simultaneously. A question then arises: is there a unified picture to describe thermoelectric properties? We note that σ, S, and κ e can be expressed in terms of integrals of a single function, the transport distribution function (TDF), which consists of carrier velocity, relaxation time, and density of states (DOS). Therefore, the essential step to search for good thermoelectrics is to find the optimal TDF, which gives high thermoelectric figure of merit. Although the transport distribution function in the delta function form has been shown to give the best thermoelectric performance mathematically [3], the shape effect of TDF in other mathematical form on thermoelectric properties has not been well studied. Moreover, in real materials, the TDF is usually having complex and broadened shape, rather than being an ideal delta function with a sharp peak. The universality of the shape effect of the TDF on thermoelectric performance becomes more and more important and emergent for the search of or the optimal design for efficient thermoelectric materials.
In this work, a systematical investigation of the shape effects of TDF in a variety of typical functional forms on thermoelectric properties of materials is carried out, in order to provide a unified picture with emphasis on asymmetrical TDF, characterized by skewness, which describes the universal features of thermoelectric power, power factor, and ZT, so as to enrich our understanding of electronic origin for improving thermoelectric efficiency. The basic theory of thermoelectricity based on the Boltzmann transport equa-tion is briefly introduced first, and the effects of TDF in terms of location parameter and scale parameter on thermoelectric performance are then presented. After that, detailed discussions are provided focusing on the asymmetric feature of TDF, characterized by skewness, as an important key to describe the thermoelectric transport properties. We have proposed a novel skewness, named as thermoelectric skewness, constructed by the symmetric and asymmetric functions of TDF, and compared our thermoelectric skewness with conventional skewness (third moment), then presented its applications combined with ab initio results and experiments.

Theory of Thermoelectricity Based on the Boltzmann Transport Equation
By solving the Boltzmann transport equation (BTE), the three properties σ, S, and κ e mentioned above can be written as integrals of a single function, the TDF [3,52]: where σðεÞ is the TDF, f μ ðT, εÞ is the Fermi-Dirac distribution function, μ is the chemical potential, and T is the temperature. The derivative of f μ ðT, εÞ is expressed as By defining the dimensionless variable x = ðε − μÞ/ðk B TÞ, constants σ 0 = e 2 /ðħa 0 Þ~46,000 (Ω cm) −1 and S 0 = k B /ẽ 87 μV K −1 , putting other material-dependent parameters into one single variable, λ = κ l /½Tσ 0 ðk B /eÞ 2 = κ l /½Tσ 0 ðS 0 Þ 2 , and using a dimensionless TDF: One can express the transport coefficients in terms of dimensionless variable and dimensionless integrals.
The reduced electrical conductivity σ r can be written as The reduced thermoelectric power (Seebeck coefficient) The reduced power factor P r is and the reduced electronic thermal conductivity κ e,r can be expressed as where the dimensionless integrals J n can be written as and the distribution density function (DDF), D n ðxÞ, is defined by The shape of D n ðxÞ is illustrated in Figure 1. It is clearly seen that D 0 ðxÞ has a maximum at x = 0 and then decreases when x is larger or smaller than 0. It drops quickly to near zero when x~±5, i.e., ðε − μÞ~±5k B T.
Finally, the figure of merit can be obtained subsequently: If we further express λ in terms of the reduced temperature T r = T/T 0 with room temperature T 0 = 300 K, the reduced lattice thermal conductivity κ l,r = κ l /κ 0 with κ 0 = 1 Wm -1 K -1 , and the dimensionless constant β 0 = κ 0 /½T 0 σ 0 ðS 0 Þ 2 = 0:0976, one can see that λ is actually a dimensionless quantity, too: For T = room temperature (RT), T r = 1, we have λ = β 0 κ l,r = 0:0976κ l,r . The expression of ZT in Equation (13) is different from the conventional expression like Equation (1). The advantage of expressing ZT in the form of Equation (13) is that every component in the equation remains dimensionless through using reduced quantity. Thus, the formula of ZT is in a neat and simple format. In this way, it is obvious that ZT is only determined by the material-dependent parameter λ (phonon part) and the dimensionless TDF gðxÞ (electron part). Because the reduced lattice thermal conductivity κ l,r (hence λ) is mainly determined by lattice contribution, the electronic mechanism for improving ZT is to optimize TDF gðxÞ. The shape of gðxÞ is therefore an important key factor for optimizing thermoelectric properties. In the following, we will examine the shape effects of TDF gðxÞ in several typical distribution forms on thermoelectric performance. We will try to link the form of TDF as closely with real materials as possible in order to enhance the applicability of the theory.

The Effects of Location and Scale
Parameters of TDF on Thermoelectric Performance The location parameter of a probability distribution in statistics is a scalar-or vector-valued parameter x 0 , which determines the shift of the distribution. For TDF, the location parameter is associated with chemical potential μ, because the energy in TDF is relative to chemical potential (ε − μ). In the dimensionless expression of TDF, the location parameter is expressed in the variance of x. The scale parameter is a special kind of numerical parameter of a parametric family of probability distributions, in probability theory and statistics. The larger the scale parameter, the more spread out the distribution. The inverse scale parameter, or the so-called rate parameter, which is the reciprocal of the scale parameter, is often used in some families of distributions. The location parameter and scale parameter are also considered together, forming a location-scale family, which is a family of probability distributions parametrized by a location parameter and a nonnegative scale parameter. The thermoelectric performance is therefore associated with the location parameter (peak positions or modes) and scale parameter of TDF. The simplest extremely case is that the TDF is only depending on scale parameter, a constant g 0 , namely, gðxÞ = g 0 , but there is no location parameter. For a constant TDF with only scale parameter, the electrical conductivity is proportional to the scale parameter, i.e., the constant g 0 , but the thermoelectric components such as S r and P r are all zero, and thus, no thermoelectric effect can be observed (see Supplemental materials (section S1) for more details).
A typical but extreme example for examining effects of location-scale family of TDF on thermoelectric properties is a single peak with Dirac delta function form. For the normalized case of TDF in the function form of Dirac delta function, the thermoelectric transport functions are determined by two constants, namely, the location parameter x 0 and the scale parameter, as had been shown in the earlier work of Mahan and Sofo [3], as well as in the recent work of Fan et al. [53]. The TDF in the single Dirac delta function form can be expressed as gðxÞ = g 0 δðx − x 0 Þ, where g 0 is constant. This TDF is actually a delta shape peak located at the x = x 0 position away from the Fermi energy level. The reduced electronic thermal conductivity is 0, independent on location-scale family (see Supplemental materials (section S2) for more details). The reduced thermoelectric power is solely depending on the location parameter, x 0 , but it is independent on the scale parameter, g 0 , while the rest of the thermoelectric components, such as the reduced electric conductivity, the reduced power factor, and the ZT are up to the scale parameter g 0 and the 3 Research location parameter x 0 , as well as the function of location parameter, D 0 ðx 0 Þ.
For relating the delta-function-like TDF to real materials, the rare-earth metal compounds are often referred to due to their very narrow f orbitals [3]. However, here, it is important to point out that cautions should be taken when linking to rare-earth metal compounds. As mentioned at the beginning of this article, the Boltzmann transport equation is used to describe the thermoelectric transport, which indicates that the systems with diffusive transport are considered, while for rare-earth metal compounds, if the narrow f orbitals are involved in thermoelectric transport, the localized orbital with possible strongly correlated features should be taken into account. In other words, to accurately describe electron strongly correlated system, the Kubo theory [54] should be used.
In real materials with itinerant electrons, the TDF usually is not as localized as delta function; it may have finite height with finite spread. Therefore, the TDF with typical general location-scale family can be described with Gaussian function form, The location parameter x 0 defines the peak position of the Gaussian function, and the scale parameter b characterizes the width of the peak. Fan et al. [53] had discussed how to optimize the ZT with TDF being in the Gaussian function form with normalized case or bounded condition.
Based on the investigation of normalized Gaussian TDF, the good thermoelectrics should have the TDF with the characteristics of two factors (see Supplemental materials (section S3) for more details). One is the optimal location parameter, which should be about x 0 = 2:4. This equally means that the peak of TDF related to Fermi energy or chemical potential should be about 2:4k B T (~60 meV at room temperature). The other factor is the scale parameter, for which, the smaller the b, the better the thermoelectric performance. This is indeed the optimization strategy for achieving high thermoelectric performance with the condition of assuming TDF is normalized. The normalized Gaussian TDF can provide some guidance for searching better thermoelectrics. However, it is not fully valid because in the real materials, the DOS is not normalized; rather, the DOS is bounded by the number of electrons in the orbitals. With the bounded condition, we have shown previously that the optimal TDF is a rectangular-shape distribution [53]. For both normalized and bounded cases, the location parameter x 0 is still valid, but the scale parameter is depending on the normalizing condition. Obviously, we need to find other factor that can characterize the thermoelectric performance in a universal manner and cover all the aforementioned aspects.
Although we try to make Gaussian TDF as general as possible, there is still something missing. The location parameter and scale parameter are describing the center and spread of a distribution. The characterization of the shape (or form) of distribution, which requires higher order of moments, is missing.
We note that, in real materials, the TDF above Fermi level (or chemical potential) can be quite different from the TDF below the Fermi level. This is actually quite a popular case for materials ranging from metal to semiconductor. In fact, the optimal TDFs such as normalized Dirac delta function [3,53] or bounded rectangle TDF [53] are asymmetric about the chemical potential. Therefore, we propose in this work that the asymmetry of the TDF may be more general and possible

Research
to be the universal characterization for thermoelectric properties. The asymmetry of the TDF is associated with the shape (or form) of distribution.

Skewness as a Key to Enhance Thermoelectric Performance
The shape of a probability distribution can be described by socalled shape parameter or form parameter, which is a kind of numerical parameter of a parametric family of probability distributions [55]. It is often measured quantitatively by the method of moments. The skewness (3 rd moment) or kurtosis (4 th moment) is often used to describe the shape. Especially, the skewness is a measure of the asymmetry of the probability distribution of a real-valued random variable about its mean. The typical function forms of TDF, including constant TDF [3], delta function TDF [3,51], Gaussian TDF [53], and rectangle TDF [53], are all symmetric and having zero skewness about their mean. To discuss TDF in a more general manner, we need to choose typical distribution function with adjustable higher order moments, especially the skewness, to describe the effects of asymmetry of the TDF on the thermoelectric properties. We found that the beta distribution function is well suited for this purpose.
The TDF using beta distribution can be expressed as where Bðα, βÞ = ðΓðαÞΓðβÞÞ/ðΓðα + βÞÞ, which is a normalization constant to ensure that the total probability is 1. Here, ΓðzÞ is the Gamma function. The beta distribution is more general than other aforementioned TDFs, and it may take a wide variety of different shapes depending on the values of the two parameters α and β within the range of [0,1], as shown in Figure 2. It can have symmetric shapes when α = β (see Figure 2(a)), taking U-shaped distribution (α = β < 1), 2-point Bernoulli distribution with equal probability 1/2 at each Dirac delta function end x = 0 and x = 1 (α = β ⟶ 0), the arcsine distribution (α = β = 1/2), uniform distribution (α = β = 1), the semielliptic distribution (α = β = 3/2), the parabolic distribution (α = β = 2), the bell-shaped distribution (α = β > 2), and the 1-point degenerate distribution with a Dirac delta function spike at the midpoint x = 0:5 with probability 1 and zero probability everywhere else (α = β ⟶ ∞). For the cases of α ≠ β (see Figure 2(b)), the beta density function is skewed. It can take a U-shaped distribution or J-shaped distribution with positive shew for α < β and negative shew for α > β. The beta distribution has wide applications in statistics, physics, project management, and so on. We have successfully applied a beta-like distribution as the statistic weight to account for the possibility of the l short region ordered structure taking place in the semiconductor alloy A x B ð1−xÞ . The band energy and the band offset as well as the formation energy for disorder alloys have been calculated using cluster expansion with above beta-like distribution [56][57][58][59][60].
The thermoelectric properties are thus can be obtained using the integrations listed in equations (7)- (13), in terms of parameters α and β. The skewness of the beta distribution can be expressed with third standardized moment, where μ m is the mean, SD is the standard deviation, and E is the expectation operator. The skewness can be positive or negative.
As we discussed before, the thermoelectric TDF has limited effective energy range about x~±5, namely, ðε − μÞ ±5k B T (see for example, Figure 1), the finite range [0,1] of beta function can be scaled and shifted to be within the 5 Research effective energy range for thermoelectric TDF. For convenient illustration, we have scaled beta TDF by 2.4, so that the case of delta function located at x = 0 or x = 1 can be reproduced. The chemical potential μ/ðk B TÞ is set to be the middle of interval, 0.5.
Firstly, we investigate the thermoelectric properties for symmetric TDF. For the symmetric beta distribution, the parameters α = β; thus, the mean = 1/2, and the skewness equals to zero. The calculated thermoelectric power (Seebeck coefficient), power factor, and ZT are found to be all zero, closely related with skewness. This is very interesting and has been shown for the first time that, for any form of TDF, if it is symmetric about its mean, the thermoelectric power (Seebeck coefficient), power factor, and ZT are zero for the case of the chemical potential μ/ðk B TÞ being set to the mean of the distribution. This important conclusion can be verified through the first order of distribution density function (DDF), D 1 ðxÞ (as shown in Figure 1), and then, x] will result in zero value. Any thermoelectric quantities that are proportional to J 1 will go to zero too. Then, we consider the cases of TDF with asymmetric shapes for α ≠ β. We start from the asymmetric linear beta TDF with slope = 2, namely, gðx, α = 2, β = 1Þ = 2x that is away from uniform distribution gðx, α = 1, β = 1Þ = 1. The reduced electrical conductivity σ r ð2, 1Þ is about 0.224 (the chemical potential μ/ðk B TÞ = 0:5), which is almost the same as σ r ð1, 1Þ, but the reduced thermoelectric power (Seebeck coefficient) increases from zero to S r ð2, 1Þ~0.082; the power factor and ZT are enhanced from 0 to 0.030 and 0.179, respectively. Similarly, a significant magnitude of skewness is found to increase from 0 to 0.566, indicating some relationship between the thermoelectric power and the skewness of the TDF. A negative skew -0.566 is found for gðx, α = 2, β = 1Þ, and positive skew is 0.566 for gðx, α = 1, β = 2Þ. The thermoelectric power (Seebeck coefficient) has the same sign to skewness.
We further explore the relationship between skewness and thermoelectric power and power factor as well as ZT for the TDF similar to parabolic distribution (gðx, α = 2, β = 2Þ) but with some degree of skewness (gðx, α = 2, βÞ). Figure 3 plots the thermoelectric properties together with the skewness of TDF in the form of gðx, α = 2, β ≥ 2Þ. The absolute value of Seebeck coefficient increases monotonically with parameter β. It can be clearly seen that the skewness of TDF is also increasing with parameter β monotonically, very similar with the trend of Seebeck coefficient. The power factor and ZT are all positive due to the square of the thermoelectric power term in the formula, and therefore, we compare with the same order of the skew, that is, γ 2 3 . As shown in Figure 3(c), the power factor can be fitted closely with square of skew with a factor of 0.151, namely, PF~0:151γ 2 3 . The variance of square of skew can also well capture the trend of ZT as a function of β ( Figure 3(d)). The close relationship between ZT and square of skew (γ 2 3 ) for beta TDF with more general cases (for any value of α and β) can be observed in Figure 4. This result is striking and surprising in twofolds. Firstly, this work presented for the first time that the asymmetry of TDF, characterized by skewness, is a universal parameter to describe the trend of thermoelectric properties for the shape of TDF. Secondly, we have shown for the first time that the thermoelectric transport properties including thermoelectric figure of merit, which are based on typical quantum statistics (Fermi-Dirac distributions), can be well described by typical shape parameter (skewness) for classical statistics.
Our results indicate that it is possible to choose suitable parameters from classical statistics to capture some important features of quantum statistics; thus, the quantum statistics and classical statistics can be linked or bridged by some typical parameters that characterizing the shape of distribution. This may suggest that a large amount of information or tools in classical statistics or probability theory may be possible to be applied for quantum statistics and used for optimization of important properties (such as the thermoelectric power, the power factor, and ZT) of quantum systems.

Universal Thermoelectric Skewness:
Construction and Applications

Construction of Universal Thermoelectric Skewness.
Above, we already clearly demonstrate that the skewness of TDF, characterization of asymmetry, can be closely related to thermoelectric power (Seebeck coefficient) and other thermoelectric properties. One remaining question is how to apply the concept of skewness in material science or, more specifically, thermoelectric materials. In case a full TDF is available, one may try to fit TDF with beta function, and then, the conventional skewness, as shown in equation 7 Research (18), can be applied directly to analyze the thermoelectric properties. Unfortunately, in the practice, there are lots of challenges of applying conventional or standard skewness in thermoelectric materials, such as the failure of fitting of TDF by beta function, or the skewness is undefined for available TDF (for example, for the Cauchy distribution, the mean, skewness, and kurtosis are undefined). At the same time, we also notice that there are several ways to measure skewness and their robustness are often argued. Considering the special characteristics of thermoelectric TDF (e.g., it is meaningful near Fermi energy or chemical potential at finite temperature, but it will be quickly decayed when the range goes beyond ±5k B T, as shown in Figure 1), it is necessary to construct a novel skewness or asymmetric function for thermoelectric TDF that can be conveniently used by materials science community.
Taking chemical potential as reference and taking into account both the location parameter and scale parameter, we may construct symmetric and asymmetric functions as The symmetric and asymmetric features of g S ðx, x 0 , bÞ and g A ðx, x 0 , bÞ for gðx, x 0 , bÞ being typical Gaussian are shown in Figure 5 with x 0 = 1 and b = 1. Because these functions are symmetric and asymmetric by definition, they are fundamental and universal for distributions with any form of shape.
Then, we can define the thermoelectric skewness as In case only the qualitative trend of thermoelectric skewness (e.g., the indication of improving direction of thermoelectric power) is needed, then integrals can be replaced by a finite number of sum (using chosen typical x i ), Due to the decay feature of the derivative of the Fermi-Dirac distribution function, the coefficients a i have larger weight for energy sampling position x i in the range about ±5k B T around chemical potential but approach to zero for x i far away from chemical potential. This interesting feature ensures that the finite number of sum is enough to well reproduce the thermopower but with much simpler form. The simplest case of the thermoelectric skewness is the ratio of asymmetric function and symmetric function at typical energy sampling positions, respectively, SymG (x) AsymG (x) Figure 5: The symmetric and asymmetric features of g S ðx, x 0 , bÞ and g A ðx, x 0 , bÞ. The gðx, x 0 , bÞ in typical Gaussian form (x 0 = 1 and b = 1), labeled as "Gaussian(x)", is also shown for comparison. The symmetric function is labeled as "SymG(x)", and the asymmetric function is labeled as "AsymG(x)", respectively.

Research
where the first approximated parameters x S1 and x A1 are the values at which D 0 ðxÞ and D 1 ðxÞ reach their maximum, respectively, namely, x S1 = 0 and x A1 = 1:543405. Surprisingly, this simplified thermoelectric skewness, γ TE Approx,1 ðx 0 , bÞ, will be shown later to have similar characteristics of third moment of skewness with mean replaced by chemical potential. Therefore, the proposed thermoelectric skewness in this work has many advantages. It is not only having characteristics of conventional skewness but also convenient to be applied in theoretical and experimental TDF for thermoelectric materials.

The Universal Thermoelectric Skewness vs. Conventional
Skewness. Here, we will compare the universal thermoelectric skewness and the conventional skewness with third moment as well as the thermoelectric power, to further illustrate the unique characteristics of the universal thermoelectric skewness. At the same time, this is also the showcase for application of universal thermoelectric skewness for prediction of thermoelectric power.
Suppose we have TDF in shifted Gaussian function form with typical scale parameter b = 2 but adjustable location parameter x 0 , we can then calculate the thermoelectric power and the skewness. For the Fisher moment coefficient of skewness (equation (18)), if the skewness is defined about the mean of distribution, then the skewness of Gaussian TDF is zero. However, if we take any value (e.g., the chemical potential) that is away from the mean as the reference point, then the distribution is no longer symmetric about the reference point. The asymmetry and skewness for shifted distribution can therefore be defined accordingly. Take the chemical potential as the reference point; we can then rewrite the conventional skewness with third standardized moment as follows: where SD is the standard deviation and E is the expectation operator. The skewness γ 3,0 is the modified version of conventional skewness. For shifted Gaussian TDF, the symmetric and asymmetric functions (Equations (19) and (18)) can be also written as The thermoelectric skewness can be obtained accordingly, It is interesting to see that both conventional skewness γ 3,0 and the thermoelectric skewness γ TE Approx,1 ðx 0 , 2Þ are in good agreement with the thermoelectric power, as shown in Figure 6. Moreover, the thermoelectric skewness not only perfectly captures the main feature of conventional skewness γ 3,0 (namely, the curvature) but also is closer to thermoelectric power. This comparison confirms the unique feature of our proposed thermoelectric skewness and its special role

Application of Universal Thermoelectric Skewness.
In this section, we will show how to apply the aforementioned proposed thermoelectric skewness of TDF (namely, asymmetricity of TDF) to interpret the thermoelectric components obtained by either first principle calculations or experimental measurements. The chosen materials cover from amorphous metal alloys (e.g., a-Al 84 Fe 16 ), bulk quasicrystals (e.g., i-AlCuFe) from experiments and quasicrystal with icosahedral (1/1) approximant from ab initio calculations. These examples represent a wide range of materials and two typical kinds of characterization of thermoelectric materials, namely, thermoelectric measurements and first principle calculations. Moreover, quasicrystals are also a new class of materials with ordered but not periodic structure, having unique surfaces and properties, as well as wide applications [61][62][63][64][65][66][67][68][69][70]. Here, we will use these typical examples to illustrate the modeling of thermoelectric skewness and make comparison with ab initio calculations and experiments.
As shown by Landauro and Solbrig [63], the spectral resistivity of icosahedral bulk quasicrystals i-Al 62  10 Research can be modeled by Lorentzians. Taking self-consistent Fermi level of the perfect approximant as reference (ε sc F = 0) and using the same symbols as above cases, the Lorentzians can be expressed as where the location parameter x 0 defines the Lorentzian peak position (with reference to ε sc F ), and the height is 1/ðπbÞ. Then, the model spectral resistivity can be expressed as a combination of two Lorentzians [63] (as shown in Figure 7 The location parameters x 0i and scaling parameters b i as well as mixing factor α and unit dependent fitting parameter A are obtained by fitting to either ab initio calculations or experiments. It is found that single Lorentzian is sufficient for amorphous a-Al 84 Fe 16 (see insert figure of Figure 7(a)), while for icosahedral bulk quasicrystals i-Al 62 Cu 25.5 Fe 12.5 and hypothetical i-AlCuFe (1/1) approximant, two Lorentzians are necessary (spectra are shown in Figures 7(a) and 7(b), and detailed parameters [63] are presented in Table S1).
The model spectral resistivity (Figure 7(a)) is related to spectral conductivity (Figure 7 In order to directly compare our thermoelectric skewness and related analysis with ab initio results and experimental measurement, in this section, we will directly utilize thermoelectric transport coefficients from Equations (2)-(5) for this purpose. Caution should be made when utilizing thermoelectric skewness (Equations (21) and (22)), and the parameters x A1 and x S1 should be scaled by k B T, in case the dimensionless quantities are no longer used in thermoelectric transport coefficients. Therefore, the parameters x A1 ðTÞ and x S1 ðTÞ are becoming temperature dependent.
The thermoelectric powers of quasicrystals are then calculated by both Boltzmann transport equation (Boltz, Equations (2)-(5)) and those from thermoelectric skewness (Equation (22)) as a function of chemical potential and compared in Figure 7(c). It is clearly seen that the thermoelectric skewness is able to capture the trend of thermoelectric power perfectly with respect to the variance of chemical potential. Similar excellent agreements are also found for the cases of amorphous a-Al 84 Fe 16 and hypothetical i-AlCuFe (1/1) approximant. Especially, the locations of chemical potential of bulk quasicrystal that leads to zero thermoelectric power or maximum value, as well as the sign-changed behavior of thermoelectric power predicted by BTE, can be well reproduced by our thermoelectric skewness.
Finally, we compare our calculated temperaturedependent thermoelectric power by BTE and thermoelectric skewness with experiments [64], as shown in Figure 7(d). In this calculation, the temperature-dependent chemical potential from ref. [63] is adapted, where ε F is Fermi energy at 0 K for i-Al 62 Cu 25.5 Fe 12.5 (labeled as "A") and i-Al 63 Cu 25 Fe 12 (labeled as "B"). Using the parameters ε F ðAÞ = 0:223 eV, ε F ðBÞ = 0:242 eV, and ξ = 0:10 × 10 −6 eV/K 2 , the agreement between theory and experiment is reasonably well. It further confirms that the proposed thermoelectric skewness is not just an alternative definition of skewness for statistics, but more importantly, it can well capture the key features of thermoelectric power, although only simplified analytical expression is used. It is worth to note that although amorphous metal alloys and quasicrystals are used for examples in this work, in principle, our proposed thermoelectric skewness of TDF can be applied for any thermoelectric materials if full or part of TDF is available. Therefore, we consider that our thermoelectric skewness is universal.
Our results are also very useful in practical searching for better thermoelectric materials. For example, it has been shown that the band engineering through the use of the thallium impurity levels in lead telluride (PbTe) enhances the thermoelectric performance [71]. Based on our findings in this work, the enhancement of thermoelectric performance found in [71] can be well understood by comparing DOS above Fermi energy and that below Fermi energy. PbTe doped with Tl introduce an additional DOS peak above the Fermi energy, by using our asymmetric function or thermoelectric skewness, with TDF being replaced by experimental DOS, the asymmetry or skewness of DOS is enhanced, and thus, the thermoelectric performance will be increased accordingly.
In general, our work suggests that, to search for better thermoelectric materials, one should first check the asymmetric feature of TDF. If full TDF is not available, it will be also useful to check the asymmetry of the components of TDF, such as DOS, relaxation time, and carrier velocity for the first screening, if they are available, either by measurements or computations.

Conclusion
In conclusion, we have systematically investigated the relationship between the parameters of TDFs (including location parameter, scale parameter, and shape parameter) and thermoelectric properties. We found that the asymmetry of TDF, characterized by skewness, is a more general and universal parameter to capture the trend of thermoelectric power, power factor as well as ZT. By taking TDF being beta distribution as a showcase, it is found that the Seebeck coefficient varies closely with the first order of skewness, while the power factor and ZT are well fitted with the second order of skewness (namely, the square of skewness). The finding of close relationship between the shape of TDF and thermoelectric performance is striking and surprising. Firstly, for any form of TDF, if it is symmetric about its mean, then the thermoelectric power (Seebeck 11 Research coefficient), power factor, and ZT are zero for the case of the chemical potential μ/ðk B TÞ being set to the mean of the distribution. Secondly, the asymmetry of TDF, characterized by skewness, is a universal parameter to describe the trend of thermoelectric properties for the shape of TDF. Thirdly, the thermoelectric transport properties including thermoelectric figure of merit, which are based on typical quantum statistics (Fermi-Dirac distributions), can be well described by typical shape parameter (skewness) for classical statistics.
In many cases, the skewness is undefined (for example, the Cauchy distribution or Lorentzian distribution) or difficult for use in materials science. It is necessary to construct a new skewness that can be adapted for materials science. We then proposed a novel skewness, so-called thermoelectric skewness, constructed from symmetric and asymmetric TDF functions. It is demonstrated that the proposed thermoelectric skewness is not only to perfectly capture the main feature of conventional skewness but also to predict accurately the thermoelectric power, by comparing with ab initio calculations and experiments. This comparison confirms the unique feature of our proposed thermoelectric skewness and its special role of connection between the statistics of TDF and thermoelectric properties of materials.
Our results indicate that it is possible to choose suitable parameters from classical statistics to capture some important features of quantum statistics; thus, the quantum statistics and classical statistics can be linked or bridged by some typical parameters that characterize the shape of distribution. This may suggest that a large amount of information or tools in classical statistics or probability theory may be applied for quantum statistics and used for optimization of important properties (such as the thermoelectric power, the power factor, and ZT) of quantum systems. Our results are also very useful in practical searching for better thermoelectric materials.

Data Availability
The data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest
There are no conflicts to declare.