Shifting a curve to another curve horizontally

This question has a few problems that needed some fiddling to resolve. This is not the ideal solution I'm sure, but it provided a close enough value to what was expected by the manual method (around 1.5 and 1.6).

The first roadblock is that when you shift the X values, you don't get matching y values, so calculating the residual can be tricky. I brute-forced my way through this problem by creating a huge new x array with 1000 points, then interpolated the original 2 y values on this new x array (this will come later). Therefore, when calculating the residual between the two curves, the x values will be off, but not by much.

reference_y = y1
to_shift_y = y2
expanded_x = np.logspace(np.log10(x[0]), np.log10(x[-1]), num=1000)
expanded_y_reference = np.interp(expanded_x, x, reference_y)
expanded_y_to_shift = np.interp(expanded_x, x, to_shift_y)

Then, when you shift x by some constant, you'll get two regions where there won't be equivalent x values.

original x: -------------------------------xxxx
shifted x:  xxxx-------------------------------

I created a new x array with the shift parameter, hor_shift set so some value greater than 1. Then, I found the indices where the original and shifted stop matching.

start = np.argmax(expanded_x >= expanded_x_shifted[0])
end = np.argmin(expanded_x_shifted <= expanded_x[-1])

Since these arrays are [False, False, True, True ...] and [True, True, ..., True, False, False], argmax and argmin will return the first instance where you have a different value.

Now, we have to slice our original and shifted x arrays so they have the same size, and values in common, and the same with the expanded y arrays. Pardon the long names, it's just so I don't get confused.

expanded_x_original_in_common_with_shifted = expanded_x[start:]
expanded_x_shifted_in_common_with_original = expanded_x_shifted[:end]
sliced_expanded_y_reference = expanded_y_reference[start:]
sliced_expanded_y_to_shift = expanded_y_to_shift[:end]

And last, and most importantly, we can calculate a distance between the two curves, assuming the x values are aligned.

residual = ((sliced_expanded_y_reference - sliced_expanded_y_to_shift) ** 2).sum()

By minimizing this, we can get the ideal shift.

We can compare our curves. Here, I used two values for the shift, 1.3 and 1.56, to illustrate good and bad shift values (these were found by testing different values). The vertical lines show the region in common.

Now, we can transform this process into a function and use some minimization method to find the ideal shift value. Here's what I got.

from lmfit import Parameters, minimize
par = Parameters()
# If the shift parameter is 1, you get an error
par.add('shift', value=1.1, min=1)

def min_function(par, x, reference_y, to_shift_y):
    hor_shift = par['shift'].value
    # print(hor_shift)  # <- in case you want to follow the process
    expanded_x = np.logspace(np.log10(x[0]), np.log10(x[-1]), num=1000)
    expanded_x_shifted = expanded_x * hor_shift
    start = np.argmax(expanded_x >= expanded_x_shifted[0])
    end = np.argmin(expanded_x_shifted <= expanded_x[-1])
    expanded_x_original_in_common_with_shifted = expanded_x[start:]
    expanded_x_shifted_in_common_with_original = expanded_x_shifted[:end]
    expanded_y_reference = np.interp(expanded_x, x, reference_y)
    expanded_y_to_shift = np.interp(expanded_x, x, to_shift_y)

    sliced_expanded_y_reference = expanded_y_reference[start:]
    sliced_expanded_y_to_shift = expanded_y_to_shift[:end]
    
    
    residual = ((sliced_expanded_y_reference - sliced_expanded_y_to_shift) ** 2).sum()
    return residual

minimize(min_function, par, method='nelder', args=(x, reference_y, to_shift_y))

This results in an ideal shift parameter of 1.555, confirming the initial guess. Note that you have to change the residual expression to (sliced_expanded_y_reference - sliced_expanded_y_to_shift) if you want your chisquared to match the one in the graphs.

The notations are changed in order to make more clear the matrix equations below :

y(x)=y1(x)

z(x)=y2(x)

A translation of value=c on the x logarithmic scale is equivalent to an expansion of value=b on the x linear scale because log(x)+c=log(b x) with c=log(b).

The inverse function x=f(y) have to approximately coincide with bx=f(z). So we consider the sum of residuals [f(y)-x]^2+[f(z)-b x]^2. This leads to the regression calculus below. The function f(y) is approximated with a polynomial of degree m.

With the given data the shape of the curve of f(y) is rather smooth. This suggest that a low degree m might be sufficient.

For example with m=2 the result is :

The black curve is the horizontaly translated blue curve from c=0.180 on the logarithmic scale.

Of course one can use a polynomial of higher degree. For example with m=3 we get b=1.536 and c=0.186

This numerical example is a favourable case because the curve x=f(y) has a simple shape. In case of more complicated shapes probably a bigger value of m should be necessary, with the risk of unreliable regression calculus.