发现两个向量之间的最佳/秤/移(finding the best/ scale/shift betw

我有两个矢量表示的函数f（x）和另一个矢量f（一个X + B），即F（X）的缩放和移位版本。 我想找到最好的规模和转移的因素。

*最好 - 通过最小二乘误差，最大似然等手段

有任何想法吗？

例如：

f1 = [0;0.450541598502498;0.0838213779969326;0.228976968716819;0.91333736150167;0.152378018969223;0.825816977489547;0.538342435260057;0.996134716626885;0.0781755287531837;0.442678269775446;0];
f2 = [-0.029171964726699;-0.0278570165494982;0.0331454732535324;0.187656956432487;0.358856370923984;0.449974662483267;0.391341738643094;0.244800719791534;0.111797007617227;0.0721767235173722;0.0854437239807415;0.143888234591602;0.251750993723227;0.478953530572365;0.748209818420035;0.908044924557262;0.811960826711455;0.512568916956487;0.22669198638799;0.168136111568694;0.365578085161896;0.644996661336714;0.823562159983554;0.792812945867018;0.656803251999341;0.545799498053254;0.587013303815021;0.777464637372241;0.962722388208354;0.980537136457874;0.734416947254272;0.375435649393553;0.106489547770962;0.0892376361668696;0.242467741982851;0.40610516900965;0.427497319032133;0.301874099075184;0.128396341665384;0.00246347624097456;-0.0322120242872125]

*注意，F（X）可能是不可逆的...

谢谢，

辖

Answer 1:

这里有一个简单，有效的，但也许有些幼稚的做法。

首先要确保你犯了一个通用的插通过两种功能。这样，你可以在给定的数据点之间评估两种功能。我用三次样条插值，因为这似乎一般容纳您提供流畅的功能类型（并且不需要额外的工具箱）。

然后你在评估大量的点源功能（“原始”）。也使用该号码作为在一个内联函数的参数，即作为输入X ，其中

X = [a b]

（如在ax+b ）。对于任何输入X ，这个内联函数将计算

在相同的x位置处的目标函数的函数值，但随后缩放和偏移由a和b分别。
所得到的函数值之间先前计算的源功能的那些相同的平方差，和总和。

使用该内联函数在fminsearch一些初始估计（已经由自动装置在视觉上或通过获得的一个）。对于您所提供的示例中，我使用了一些随机的，这都收敛到接近最优的拟合。

所有代码上面的：

function s = findScaleOffset

    %% initialize

    f2 = [0;0.450541598502498;0.0838213779969326;0.228976968716819;0.91333736150167;0.152378018969223;0.825816977489547;0.538342435260057;0.996134716626885;0.0781755287531837;0.442678269775446;0];
    f1 = [-0.029171964726699;-0.0278570165494982;0.0331454732535324;0.187656956432487;0.358856370923984;0.449974662483267;0.391341738643094;0.244800719791534;0.111797007617227;0.0721767235173722;0.0854437239807415;0.143888234591602;0.251750993723227;0.478953530572365;0.748209818420035;0.908044924557262;0.811960826711455;0.512568916956487;0.22669198638799;0.168136111568694;0.365578085161896;0.644996661336714;0.823562159983554;0.792812945867018;0.656803251999341;0.545799498053254;0.587013303815021;0.777464637372241;0.962722388208354;0.980537136457874;0.734416947254272;0.375435649393553;0.106489547770962;0.0892376361668696;0.242467741982851;0.40610516900965;0.427497319032133;0.301874099075184;0.128396341665384;0.00246347624097456;-0.0322120242872125];

    figure(1), clf, hold on

    h(1) = subplot(2,1,1); hold on
    plot(f1);
    legend('Original')

    h(2) = subplot(2,1,2); hold on
    plot(f2);

    linkaxes(h)
    axis([0 max(length(f1),length(f2)), min(min(f1),min(f2)),max(max(f1),max(f2))])


    %% make cubic interpolators and test points

    pp1 = spline(1:numel(f1), f1);
    pp2 = spline(1:numel(f2), f2);

    maxX = max(numel(f1), numel(f2));
    N  = 100 * maxX;

    x2 = linspace(1, maxX, N);
    y1 = ppval(pp1, x2);

    %% search for parameters

    s = fminsearch(@(X) sum( (y1 - ppval(pp2,X(1)*x2+X(2))).^2 ), [0 0])

    %% plot results

    y2 = ppval( pp2, s(1)*x2+s(2));

    figure(1), hold on
    subplot(2,1,2), hold on    
    plot(x2,y2, 'r')
    legend('before', 'after')



end

结果：

s =
2.886234493867320e-001    3.734482822175923e-001

请注意，这个计算从您所产生的数据的一个逆转换。倒车的数字：

>> 1/s(1) 
ans =    
    3.464721948700991e+000   % seems pretty decent 
>> -s(2)
ans = 
    -3.734482822175923e-001  % hmmm...rather different from 7/11!

（我不知道关于您所提供的7/11值;用你给做图的结果不太准确逼近源功能的准确值...你确定的7/11？）

准确性可以通过能够提高

使用不同的优化器（ fmincon ， fminunc等）
从要求以更高的精度fminsearch通过optimset
具有在这两个多个采样点f1和f2以改善内插的质量
使用更好的初始估计

无论如何，这种做法是相当普遍的，并给出很好的结果。它也不需要工具箱。

它虽然一个主要的缺点-发现可能不是全局优化 ，例如解决方案，这种方法的结果的质量可能是您提供的初步估计相当敏感。所以，总是让一个（差）的情节，以确保最终的解决方案是正确的，或者如果你有大量的这样的事情要做，在你决定重新开始与不同的优化计算某种品质因数初步估计。

当然，很可能使用傅立叶+ Mellin变换的结果（由下面chaohuang所建议的）作为初始估计该方法。这可能是矫枉过正，为您提供简单的例子，但我可以很容易想象的情况下，这可能确实是非常有用的。

Answer 2:

对于每一个f(x)采取的绝对值f(x)和其归一化，使得它可以被认为在它的载体上的概率质量函数。计算预期值E[x]和方差Var[x] 然后，我们有

  E[a x + b] = a E[x] + b
Var[a x + b] = a^2 Var[x]

使用上面的方程和的已知值E[x]和Var[x]以计算a和b 。以你的价值观f1和f2从你的榜样，下面的八度脚本执行此过程：

% Octave script
% f1, f2 are defined as given in your example
f1 = [zeros(length(f2) - length(f1), 1); f1];

save_f1 = f1; save_f2 = f2;

f1 = abs( f1 ); f2 = abs( f2 );
f1 = f1 ./ sum( f1 ); f2 = f2 ./ sum( f2 );

mean = @(x)sum(((1:length(x))' .* x));
var = @(x)sum((((1:length(x))'-mean(x)).^2) .* x);

m1 = mean(f1); m2 = mean(f2);
v1 = var(f1); v2 = var(f2)

a = sqrt( v2 / v1 ); b = m2 - a * m1;

plot( a .* (1:length( save_f1 )) + b, save_f1, ...
      1:length( save_f2 ), save_f2 );
axis([0 length( save_f1 )];

而输出

Answer 3:

为比例因子，可以通过由于傅立叶变换是不变的移位计算所述两个信号的振幅频谱的比估计它。

类似地，可以通过使用梅林变换，这是尺度不变估计移因子B。

Answer 4:

这里是一个超级简单的估算尺度方法a ，关于您的示例数据的工作原理：

a = length(f2) / length(f1)

这给了3.4167，接近3.4的规定值。如果这估计是不够好，你可以使用相关性来估计转变。

我知道这是不是正是你问什么，但它可能是依赖于数据的可接受的替代方法。

Answer 5:

无论罗迪Oldenhuis和jstarr的答案是正确的。我将我自己的答案只是总结东西，以及它们之间的连接。我搞砸了罗迪的代码一点点，结束了以下内容：

function findScaleShift
load f1f2

x0 = [length(f1)/length(f2) 0]; %initial guess, can do better
n=length(f1);
costFunc = @(z) sum((eval_f1(z,f2,n)-f1).^2);
opt.TolFun = eps; 
xopt=fminsearch(costFunc,x0,opt);
f1r=eval_f1(xopt,f2,n);
subplot(211);
plot(1:n,f1,1:n,f1r,'--','linewidth',5)
title(xopt);
subplot(212);
plot(1:n,(f1-f1r).^2);
title('squared error')
end


function y = eval_f1(x,f2,n)
t = maketform('affine',[x(1) 0 x(2); 0 1 0 ; 0 0 1]');
y=imtransform(f2',t,'cubic','xdata',[1 n ],'ydata',[1 1])';
end

这使零组的结果：这种方法是准确的，但详尽的，可能需要一些时间。另一个缺点是，它仅找到局部极小，并且可能产生假的结果，如果初始猜测（X0）远。

在另一方面，jstarr方法得出以下结果：

xopt = [ 3.49655562549115         -0.676062367063033]

这是正确的答案10％的偏差。相当快的解决方案，但不是我要求的，但还是应该注意准确。我认为，为了获得最好的效果jstarr方法应该被用作通过罗迪旨意方法的初始猜测，给人一种精确的解决方案。

辖

文章来源: finding the best/ scale/shift between two vectors