Demo 2: Two-parameter scaling
In this demo we will show how to perform finite-size scaling where two free parameters are tuned to obtain the best possible data collapse. As in Demo 1, we will use the CDW transition in the square lattice Holstein model as an example, but this time both the critical inverse temperature $\beta_c$ and the exponent $2 - \eta \equiv \gamma/\nu$ are to be optimized, i.e. we do not assume the exponent is equal to its theoretical value of $7/4$ for the Ising universality class.
Let us use the same data as in Demo 1, which as before has been loaded into an array of the form:
data_with_error = [[X_L4, Y_L4, E_L4, 4], [X_L6, Y_L6, E_L6, 6], [X_L8, Y_L8, E_L8, 8], [X_L10, Y_L10, E_L10, 10], [X_L12, Y_L12, E_L12, 12]]
We now aim to rescale the data along new axes $X_s = (\beta - \beta_c)L$ and $Y_s = S_{cdw} L^{-(2-\eta)}$. We will take $\beta_c$ to be the first free parameter we will tune to obtain the optimal data collapse, and denote it $v_1$, and we will take $(2-\eta)$ to be the second free parameter, and denote it $v_2$. That is, our scaled axes will be $X_s = (X - v_1)L$ and $Y_s = Y L^{-v_2}$. For two-parameter scaling, the user now must define two functions which define the scaled X and Y axes:
x_scaled(X, L, v1, v2) = (X .- v1) * L
y_scaled(Y, L, v1, v2) = Y * (L^(-1 * v2))
The functions defining the $X_s$ axis should take $X$, $L$, $v_1$, and $v_2$ as arguments, while the function defining the scaled $Y_s$ axis should take $Y$, $L$, $v_1$, and $v_2$ as arguments. Note that since $X$ and $Y$ are arrays, elementwise operators such as $.-$ may be necessary.
To determine the best value of $v_1$, a search is performed between an initial value $v_{1i}$ and a final value $v_{1f}$. To determine the best value of $v_2$, a search is performed between an initial value $v_{2i}$ and a final value $v_{2f}$. The finite-size scaling method works in the same way as for one-parameter scaling, except now the method finds the pair of $v_1$ and $v_2$ values which yield the minimal overall residual. As in Demo 1, we will use inverse variance weights to perform weighted least-squares.
fit_weights = [1 ./ (E_L4.^2), 1 ./ (E_L6.^2), 1 ./ (E_L8.^2), 1 ./ (E_L10.^2), 1 ./ (E_L12.^2)]
However, note that since our scaled $Y_s$ axis depends exponentially on the tuned parameter $v_2$, it is important to "normalize" each calculated residual to obtain the best data collapse. To do this we will set norm_y=true
when calling the two-parameter finite-size scaling function fss_two_var
.
When the scaled vertical axis $Y_s$ has an explicit dependence on the tuned parameter $v_1$ and/or $v_2$, the magnitude of $Y_s$ values can vary drastically during the parameter sweep. In determining the relative quality of fit of a polynomial to the scaled data, the absolute magnitude of the sum of squared residuals therefore may not be an appropriate choice.
In this case, each individual fit residual should be divided by the magnitude $Y_s$ of the corresponding data point. Then, the sum of the squares of these "normalized" residuals will provide a true measure of the relative quality of the polynomial fitting.
To do this, both fss_one_var
and fss_two_var
take an optional boolean argument norm_y
. Set norm_y=true
to normalize the fit residuals.
The next step is to call the function fss_two_var
to perform the finite-size scaling.
julia> scaled_data, residuals, min_res, best_v1, best_v2 = fss_two_var(data=data_with_error, xs=x_scaled, ys=y_scaled, v1i=5.0, v1f=7.0, n1=100, v2i=1.0, v2f=2.0, n2=100, p=4, weights=fit_weights, norm_y=true);
where data
is the single array of data defined previously, and xs
and ys
are the functions previously defined for the scaled axes. Here v1i
and v1f
are the start and end points of the parameter search for $v_1$ (in this example $\beta_c$), where n1
is the number of values of $v_1$ in this range to check. Similarly, v2i
and v2f
are the start and end points of the parameter search for $v_2$ (in this example $(2-\eta)$), where n2
is the number of values of $v_2$ in this range to check. The integer degree $p$ of the polynomial must also be specified, typically $4 \leq p \leq 8$ is sufficient.
The function fss_two_var
returns five variables: an array scaled_data_array
where each element is an array of $[X_s, Y_s, E_s, L]$ data for a given lattice size; a two-dimensional array residuals
of dimension (n2, n1)
which stores the sum of squared residuals for each pair of $(v_2, v_1)$ values checked; a scalar min_res
which is the minimum value of the array residuals
; and scalars best_v1
and best_v2
which are the values of $v_1$ and $v_2$ which gave the smallest overall residual. By default it will also print out the values of best_v1
, best_v2
, and min_res
:
Optimal v1 value: 6.090909090909091
Optimal v2 value: 1.6767676767676767
Smallest residual: 1.158515739950633
At this point, you can again call the function plot_data
, passing in the scaled_data_array
returned by the fss_two_var
function. This will produce a plot of the optimal data collapse, i.e. the scaled data with $v_1$ set to best_v1
and $v_2$ set to best_v2
:
julia> plot_data(scaled_data)
Contour plot for two-parameter scaling
When performing two-parameter scaling, it is often useful to produce a contour plot showing how the quality of the data collapse varies in the $(v_1, v_2)$ plane. To make this plot, use the function plot_contour
, passing in the two-dimensional array residuals
returned by fss_two_var
. You will need to pass in the same values of v1i
, v1f
, n1
, v2i
, v2f
, and n2
which were used in fss_two_var
. You will also need to specify how many contour lines are drawn using the argument levels
. levels
can either be an integer or an array of values. If an array of values is specified, only contour lines at these specific values will be drawn on the plot. If levels
is an integer, this sets the total number of contour lines which will be drawn. The values represented by the contour lines can either be equally spaced or logarithmically spaced, which is set by the boolean argument logspace
. Often logarthmic spacing is best in order to get a reasonable number of contour lines plotted near the minima, and logspace=true
by default. A solid color fill between contour lines will be drawn if the boolean argument fill=true
, which it is by default.
julia> plot_contour(residuals, v1i=5.0, v1f=7.0, n1=100, v2i=1.0, v2f=2.0, n2=100, levels=30, fill=true, logspace=true, xlabel=L"\beta_c", ylabel=L"2 - \eta")
The plot_contour
function produces a contour plot displaying the sum of squared residuals in the $(v_1, v_2)$ plane, where a smaller value corresponds to a better quality data collapse. The minimum residual value is also indicated by a marker on the plot, showing the location of (best_v1, best_v2)
:
Note that various plot attributes such as the color scheme used, plot dimensions, marker size, marker color, marker shape, axes labels, and font sizes can be customized when calling the function plot_contour
. For example, the color
argument can be set to any color scheme supported by Plots.jl. See the docstrings of plot_contour
or the Methods page for full details. Below is a contour plot showing the same residuals
data as above, where several of these plot attributes have been modified.
julia> plot_contour(residuals, v1i=5.0, v1f=7.0, n1=100, v2i=1.0, v2f=2.0, n2=100, levels=25, fill=true, logspace=true, xlabel=L"\beta_c", ylabel=L"2 - \eta", color=:terrain, markersize=7, markershape=:star5, markercolor=:white, size=(800,300))
To use $\LaTeX$ when specifying axes labels, you can use a LaTeXString e.g. L"\beta_c"
. To use the LaTeXStrings.jl package, add it to your environment with ] add LaTeXStrings
and then import it:
julia> using LaTeXStrings