代码之家  ›  专栏  ›  技术社区  ›  Nick Fortescue

如何计算一组循环数据的平均值?

  •  130
  • Nick Fortescue  · 技术社区  · 17 年前

    我要计算一组循环数据的平均值。例如,我可能有几个罗盘读数的样本。当然,问题是如何处理包围。同样的算法可能对时钟面很有用。

    实际的问题更复杂——在一个“环绕”的球体或代数空间中,例如加法组mod n,统计学意味着什么?答案可能不是唯一的,例如,359度和1度的平均值可以是0度或180度,但统计学上0看起来更好。

    这对我来说是一个真正的编程问题,我试图让它看起来不只是一个数学问题。

    30 回复  |  直到 7 年前
        1
  •  84
  •   starblue    17 年前

    从角度计算单位向量,并取其平均值的角度。

        2
  •  51
  •   Nick Fortescue    15 年前

    这一问题在书中有详细的论述: “球体统计”,杰弗里·S·沃森,阿肯色大学讲座 《数学科学笔记》,1983年,John Wiley&Sons,Inc.如前所述 http://catless.ncl.ac.uk/Risks/7.44.html#subj4 Bruce Karsh。

    从一组角度测量值估计平均角度A的好方法 A[I] 0和lt;= i

                       sum_i_from_1_to_N sin(a[i])
    a = arctangent ---------------------------
                       sum_i_from_1_to_N cos(a[i])
    

    Starblue给出的方法在计算上是等价的,但是他的理由更清楚,而且可能在程序上更有效,而且在零情况下也能很好地工作,所以他很荣幸。

    现在对这个主题进行了更详细的探讨。 on Wikipedia 以及其他用途,如分数部分。

        3
  •  45
  •   Alnitak    17 年前

    我看到了问题——例如,如果你有一个45'角和一个315'角,那么“自然”平均值将是180',但是你想要的值实际上是0。

    我觉得Starblue有点意思。只需计算每个角度的(x,y)笛卡尔坐标,并将这些结果向量相加。最终矢量的角偏移应该是您需要的结果。

    x = y = 0
    foreach angle {
        x += cos(angle)
        y += sin(angle)
    }
    average_angle = atan2(y, x)
    

    现在我忽略了指南针的方向是从北开始,顺时针方向,而“正常”笛卡尔坐标是从0开始沿着x轴,然后逆时针方向。不管怎样,数学的计算方法应该是一样的。

        4
  •  16
  •   darron    16 年前

    对于两个角度的特殊情况:

    答案 ((A+B)360型)/2 错误的 . 对于角350和2,最近的点是356,而不是176。

    单位向量和三角解可能太贵了。

    我从一个小修补中得到的是:

    diff = ( ( a - b + 180 + 360 ) mod 360 ) - 180
    angle = (360 + b + ( diff / 2 ) ) mod 360
    
    • 0,180->90(有两个答案:该方程取a的顺时针答案)
    • 180,0->270(见上文)
    • 180,1->90.5
    • 180->90.5
    • 20, 350和5;
    • 350,20->5(以下所有示例也都正确反转)
    • 10,20->15
    • 350,2->356
    • 359,0->359.5
    • 180,180->180
        5
  •  14
  •   Adaline Simonian    8 年前

    Ackb是对的,这些基于向量的解不能被视为角度的真实平均值,它们只是单位向量对应值的平均值。然而,ACKB的建议解决方案似乎在数学上并不合理。

    以下是从最小化目标(角度[i]-Avgangle)^2(如有必要,对差异进行校正)数学推导出的一个解决方案,这使得它成为角度的真正算术平均值。

    首先,我们需要准确地看看哪些情况下角度之间的差异不同于它们的正常数对应项之间的差异。考虑角度x和y,如果y>=x-180和y<=x+180,那么我们可以直接使用差分(x-y)。否则,如果第一个条件不满足,那么我们必须在计算中使用(y+360)而不是y。相应地,如果第二个条件不满足,那么我们必须使用(y-360)而不是y。由于曲线方程,我们只在这些不等式从真变为假或从真变为假的点上最小化,因此我们可以将将整个[0360]范围划分为一组由这些点分隔的段。然后,我们只需要找到每个段的最小值,然后找到每个段最小值的最小值,即平均值。

    这是一幅图像,演示了计算角度差时出现的问题。如果x位于灰色区域,则会出现问题。

    Angle comparisons

    为了最小化一个变量,根据曲线,我们可以取我们想要最小化的导数,然后找到转折点(导数=0)。

    在这里,我们将应用最小化平方差的思想来推导常见的算术平均公式:和(a[i])/n。曲线y=和((a[i]-x)^2)可以这样最小化:

    y = sum((a[i]-x)^2)
    = sum(a[i]^2 - 2*a[i]*x + x^2)
    = sum(a[i]^2) - 2*x*sum(a[i]) + n*x^2
    
    dy\dx = -2*sum(a[i]) + 2*n*x
    
    for dy/dx = 0:
    -2*sum(a[i]) + 2*n*x = 0
    -> n*x = sum(a[i])
    -> x = sum(a[i])/n
    

    现在将其应用于具有我们调整后差异的曲线:

    b=a的子集,其中正确(角度)差a[i]-x C=a的子集,其中正确(角)差(a[i]-360)-x C= C的大小 d=a的子集,其中正确(角度)差(a[i]+360)-x d=d的大小

    y = sum((b[i]-x)^2) + sum(((c[i]-360)-b)^2) + sum(((d[i]+360)-c)^2)
    = sum(b[i]^2 - 2*b[i]*x + x^2)
      + sum((c[i]-360)^2 - 2*(c[i]-360)*x + x^2)
      + sum((d[i]+360)^2 - 2*(d[i]+360)*x + x^2)
    = sum(b[i]^2) - 2*x*sum(b[i])
      + sum((c[i]-360)^2) - 2*x*(sum(c[i]) - 360*cn)
      + sum((d[i]+360)^2) - 2*x*(sum(d[i]) + 360*dn)
      + n*x^2
    = sum(b[i]^2) + sum((c[i]-360)^2) + sum((d[i]+360)^2)
      - 2*x*(sum(b[i]) + sum(c[i]) + sum(d[i]))
      - 2*x*(360*dn - 360*cn)
      + n*x^2
    = sum(b[i]^2) + sum((c[i]-360)^2) + sum((d[i]+360)^2)
      - 2*x*sum(x[i])
      - 2*x*360*(dn - cn)
      + n*x^2
    
    dy/dx = 2*n*x - 2*sum(x[i]) - 2*360*(dn - cn)
    
    for dy/dx = 0:
    2*n*x - 2*sum(x[i]) - 2*360*(dn - cn) = 0
    n*x = sum(x[i]) + 360*(dn - cn)
    x = (sum(x[i]) + 360*(dn - cn))/n
    

    仅此一项还不足以获得最小值,但它适用于具有无边界集的正常值,因此结果肯定在集的范围内,因此是有效的。我们需要一个范围内的最小值(由段定义)。如果最小值小于段的下限,则该段的最小值必须位于下限(因为二次曲线只有1个转折点),如果最小值大于段的上限,则该段的最小值位于上限。当我们有了每个段的最小值后,我们就可以简单地找到最小值的那个(sum((b[i]-x)^2)+sum((c[i]-360)-b)^2)+sum((d[i]+360)-c)^2))。

    这是曲线的图像,它显示了在x=(a[i]+180)%360的点上它是如何变化的。有问题的数据集是65,92230320250。

    Curve

    这里是Java中算法的实现,包括一些优化,其复杂度为O(nLogn)。如果将基于比较的排序替换为非基于比较的排序(如基数排序),则可以将其简化为O(N)。

    static double varnc(double _mean, int _n, double _sumX, double _sumSqrX)
    {
        return _mean*(_n*_mean - 2*_sumX) + _sumSqrX;
    }
    //with lower correction
    static double varlc(double _mean, int _n, double _sumX, double _sumSqrX, int _nc, double _sumC)
    {
        return _mean*(_n*_mean - 2*_sumX) + _sumSqrX
                + 2*360*_sumC + _nc*(-2*360*_mean + 360*360);
    }
    //with upper correction
    static double varuc(double _mean, int _n, double _sumX, double _sumSqrX, int _nc, double _sumC)
    {
        return _mean*(_n*_mean - 2*_sumX) + _sumSqrX
                - 2*360*_sumC + _nc*(2*360*_mean + 360*360);
    }
    
    static double[] averageAngles(double[] _angles)
    {
        double sumAngles;
        double sumSqrAngles;
    
        double[] lowerAngles;
        double[] upperAngles;
    
        {
            List<Double> lowerAngles_ = new LinkedList<Double>();
            List<Double> upperAngles_ = new LinkedList<Double>();
    
            sumAngles = 0;
            sumSqrAngles = 0;
            for(double angle : _angles)
            {
                sumAngles += angle;
                sumSqrAngles += angle*angle;
                if(angle < 180)
                    lowerAngles_.add(angle);
                else if(angle > 180)
                    upperAngles_.add(angle);
            }
    
    
            Collections.sort(lowerAngles_);
            Collections.sort(upperAngles_,Collections.reverseOrder());
    
    
            lowerAngles = new double[lowerAngles_.size()];
            Iterator<Double> lowerAnglesIter = lowerAngles_.iterator();
            for(int i = 0; i < lowerAngles_.size(); i++)
                lowerAngles[i] = lowerAnglesIter.next();
    
            upperAngles = new double[upperAngles_.size()];
            Iterator<Double> upperAnglesIter = upperAngles_.iterator();
            for(int i = 0; i < upperAngles_.size(); i++)
                upperAngles[i] = upperAnglesIter.next();
        }
    
        List<Double> averageAngles = new LinkedList<Double>();
        averageAngles.add(180d);
        double variance = varnc(180,_angles.length,sumAngles,sumSqrAngles);
    
        double lowerBound = 180;
        double sumLC = 0;
        for(int i = 0; i < lowerAngles.length; i++)
        {
            //get average for a segment based on minimum
            double testAverageAngle = (sumAngles + 360*i)/_angles.length;
            //minimum is outside segment range (therefore not directly relevant)
            //since it is greater than lowerAngles[i], the minimum for the segment
            //must lie on the boundary lowerAngles[i]
            if(testAverageAngle > lowerAngles[i]+180)
                testAverageAngle = lowerAngles[i];
    
            if(testAverageAngle > lowerBound)
            {
                double testVariance = varlc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,i,sumLC);
    
                if(testVariance < variance)
                {
                    averageAngles.clear();
                    averageAngles.add(testAverageAngle);
                    variance = testVariance;
                }
                else if(testVariance == variance)
                    averageAngles.add(testAverageAngle);
            }
    
            lowerBound = lowerAngles[i];
            sumLC += lowerAngles[i];
        }
        //Test last segment
        {
            //get average for a segment based on minimum
            double testAverageAngle = (sumAngles + 360*lowerAngles.length)/_angles.length;
            //minimum is inside segment range
            //we will test average 0 (360) later
            if(testAverageAngle < 360 && testAverageAngle > lowerBound)
            {
                double testVariance = varlc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,lowerAngles.length,sumLC);
    
                if(testVariance < variance)
                {
                    averageAngles.clear();
                    averageAngles.add(testAverageAngle);
                    variance = testVariance;
                }
                else if(testVariance == variance)
                    averageAngles.add(testAverageAngle);
            }
        }
    
    
        double upperBound = 180;
        double sumUC = 0;
        for(int i = 0; i < upperAngles.length; i++)
        {
            //get average for a segment based on minimum
            double testAverageAngle = (sumAngles - 360*i)/_angles.length;
            //minimum is outside segment range (therefore not directly relevant)
            //since it is greater than lowerAngles[i], the minimum for the segment
            //must lie on the boundary lowerAngles[i]
            if(testAverageAngle < upperAngles[i]-180)
                testAverageAngle = upperAngles[i];
    
            if(testAverageAngle < upperBound)
            {
                double testVariance = varuc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,i,sumUC);
    
                if(testVariance < variance)
                {
                    averageAngles.clear();
                    averageAngles.add(testAverageAngle);
                    variance = testVariance;
                }
                else if(testVariance == variance)
                    averageAngles.add(testAverageAngle);
            }
    
            upperBound = upperAngles[i];
            sumUC += upperBound;
        }
        //Test last segment
        {
            //get average for a segment based on minimum
            double testAverageAngle = (sumAngles - 360*upperAngles.length)/_angles.length;
            //minimum is inside segment range
            //we test average 0 (360) now           
            if(testAverageAngle < 0)
                testAverageAngle = 0;
    
            if(testAverageAngle < upperBound)
            {
                double testVariance = varuc(testAverageAngle,_angles.length,sumAngles,sumSqrAngles,upperAngles.length,sumUC);
    
                if(testVariance < variance)
                {
                    averageAngles.clear();
                    averageAngles.add(testAverageAngle);
                    variance = testVariance;
                }
                else if(testVariance == variance)
                    averageAngles.add(testAverageAngle);
            }
        }
    
    
        double[] averageAngles_ = new double[averageAngles.size()];
        Iterator<Double> averageAnglesIter = averageAngles.iterator();
        for(int i = 0; i < averageAngles_.length; i++)
            averageAngles_[i] = averageAnglesIter.next();
    
    
        return averageAngles_;
    }
    

    一组角度的算术平均值可能不符合你对平均值应该是多少的直觉。例如,集合179179、0181181的算术平均值为216(和144)。您马上想到的答案可能是180,但是众所周知,算术平均值受到边值的严重影响。你还应该记住,角度不是向量,就像有时处理角度时看起来那样吸引人。

    当然,该算法也适用于所有遵循模块化算法的数量(最小调整),如一天中的时间。

    我还要强调的是,尽管这是角度的真实平均值,但与向量解不同,这并不一定意味着它是您应该使用的解,相应单位向量的平均值很可能是您实际应该使用的值。

        6
  •  6
  •   David Hanak    17 年前

    你必须定义 平均的 更准确地说。对于两个角度的具体情况,我可以考虑两种不同的情况:

    1. “真实”平均值,即(a+b)/2%360。
    2. 当保持在同一个半圆内时,在其他两个半圆之间指向“中间”的角度,例如,对于355和5,这将是0,而不是180。为此,您需要检查两个角度之间的差异是否大于180。如果是这样,在使用上述公式之前,将较小的角度增加360。

    不过,我不知道第二种方法如何推广到两个以上的角度。

        7
  •  4
  •   Community Mohan Dere    15 年前

    和所有的平均值一样,答案取决于度量的选择。对于给定的公制m,在[1,n]中,k的某些角a_k in[-pi,pi]的平均值是角a_m,它最小化平方距离d^2_m(a_m,a_k)之和。对于加权平均数,只需在总和中包括权重w_k(例如,sum_k w_k=1)。也就是说,

    a_m=arg min_x sum_k w_k d^2_m(x,a_k)

    度量的两个常见选择是Frobenius和Riemann度量。对于弗罗贝尼乌斯度量,存在一个直接公式,对应于圆统计中常用的平均轴承概念。参见“旋转组的平均值和平均值”,Maher Moakher,《暹罗矩阵分析与应用杂志》,第24卷,2002年第1期,了解详细信息。
    http://link.aip.org/link/?SJMAEL/24/1/1

    下面是GNU八度音阶3.2.4的一个函数,用于计算:

    function ma=meanangleoct(a,w,hp,ntype)
    %   ma=meanangleoct(a,w,hp,ntype) returns the average of angles a
    %   given weights w and half-period hp using norm type ntype
    %   Ref: "Means and Averaging in the Group of Rotations",
    %   Maher Moakher, SIAM Journal on Matrix Analysis and Applications,
    %   Volume 24, Issue 1, 2002.
    
    if (nargin<1) | (nargin>4), help meanangleoct, return, end 
    if isempty(a), error('no measurement angles'), end
    la=length(a); sa=size(a); 
    if prod(sa)~=la, error('a must be a vector'); end
    if (nargin<4) || isempty(ntype), ntype='F'; end
    if ~sum(ntype==['F' 'R']), error('ntype must be F or R'), end
    if (nargin<3) || isempty(hp), hp=pi; end
    if (nargin<2) || isempty(w), w=1/la+0*a; end
    lw=length(w); sw=size(w); 
    if prod(sw)~=lw, error('w must be a vector'); end
    if lw~=la, error('length of w must equal length of a'), end
    if sum(w)~=1, warning('resumming weights to unity'), w=w/sum(w); end
    
    a=a(:);     % make column vector
    w=w(:);     % make column vector
    a=mod(a+hp,2*hp)-hp;    % reduce to central period
    a=a/hp*pi;              % scale to half period pi
    z=exp(i*a); % U(1) elements
    
    % % NOTA BENE:
    % % fminbnd can get hung up near the boundaries.
    % % If that happens, shift the input angles a
    % % forward by one half period, then shift the
    % % resulting mean ma back by one half period.
    % X=fminbnd(@meritfcn,-pi,pi,[],z,w,ntype);
    
    % % seems to work better
    x0=imag(log(sum(w.*z)));
    X=fminbnd(@meritfcn,x0-pi,x0+pi,[],z,w,ntype);
    
    % X=real(X);              % truncate some roundoff
    X=mod(X+pi,2*pi)-pi;    % reduce to central period
    ma=X*hp/pi;             % scale to half period hp
    
    return
    %%%%%%
    
    function d2=meritfcn(x,z,w,ntype)
    x=exp(i*x);
    if ntype=='F'
        y=x-z;
    else % ntype=='R'
        y=log(x'*z);
    end
    d2=y'*diag(w)*y;
    return
    %%%%%%
    
    % %   test script
    % % 
    % % NOTA BENE: meanangleoct(a,[],[],'R') will equal mean(a) 
    % % when all abs(a-b) < pi/2 for some value b
    % % 
    % na=3, a=sort(mod(randn(1,na)+1,2)-1)*pi;
    % da=diff([a a(1)+2*pi]); [mda,ndx]=min(da);
    % a=circshift(a,[0 2-ndx])    % so that diff(a(2:3)) is smallest
    % A=exp(i*a), B1=expm(a(1)*[0 -1; 1 0]), 
    % B2=expm(a(2)*[0 -1; 1 0]), B3=expm(a(3)*[0 -1; 1 0]),
    % masimpl=[angle(mean(exp(i*a))) mean(a)]
    % Bsum=B1+B2+B3; BmeanF=Bsum/sqrt(det(Bsum)); 
    % % this expression for BmeanR should be correct for ordering of a above
    % BmeanR=B1*(B1'*B2*(B2'*B3)^(1/2))^(2/3);
    % mamtrx=real([[0 1]*logm(BmeanF)*[1 0]' [0 1]*logm(BmeanR)*[1 0]'])
    % manorm=[meanangleoct(a,[],[],'F') meanangleoct(a,[],[],'R')]
    % polar(a,1+0*a,'b*'), axis square, hold on
    % polar(manorm(1),1,'rs'), polar(manorm(2),1,'gd'), hold off
    
    %     Meanangleoct Version 1.0
    %     Copyright (C) 2011 Alphawave Research, robjohnson@alphawaveresearch.com
    %     Released under GNU GPLv3 -- see file COPYING for more info.
    %
    %     Meanangle is free software: you can redistribute it and/or modify
    %     it under the terms of the GNU General Public License as published by
    %     the Free Software Foundation, either version 3 of the License, or (at
    %     your option) any later version.
    %
    %     Meanangle is distributed in the hope that it will be useful, but
    %     WITHOUT ANY WARRANTY; without even the implied warranty of
    %     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
    %     General Public License for more details.
    %
    %     You should have received a copy of the GNU General Public License
    %     along with this program.  If not, see `http://www.gnu.org/licenses/'.
    
        8
  •  3
  •   Tom K    13 年前

    我想和一个没有浮点或三角学功能的微控制器分享我使用的方法。我仍然需要“平均”10个原始轴承读数,以消除变化。

    1. 检查第一个轴承的范围是270-360度还是0-90度(北部两个象限)
    2. 如果是,将此读数和所有后续读数旋转180度,将所有值保持在0<=轴承<360范围内。否则,在读数出现时进行测量。
    3. 一旦读取10个读数,计算数值平均值,假设没有包围。
    4. 如果180度旋转生效,则将计算出的平均值旋转180度,以返回“真实”轴承。

    这不理想,它可能会破裂。在这种情况下,我成功了,因为这个装置的旋转速度非常慢。我会把它放在那里,以防其他人发现自己在类似的限制下工作。

        9
  •  3
  •   DuduArbel    11 年前

    以下是完整的解决方案: (输入是以度为单位的方位数组(0-360)

    public static int getAvarageBearing(int[] arr)
    {
        double sunSin = 0;
        double sunCos = 0;
        int counter = 0;
    
        for (double bearing : arr)
        {
            bearing *= Math.PI/180;
    
            sunSin += Math.sin(bearing);
            sunCos += Math.cos(bearing);
            counter++; 
        }
    
        int avBearing = INVALID_ANGLE_VALUE;
        if (counter > 0)
        {
            double bearingInRad = Math.atan2(sunSin/counter, sunCos/counter);
            avBearing = (int) (bearingInRad*180f/Math.PI);
            if (avBearing<0)
                avBearing += 360;
        }
    
        return avBearing;
    }
    
        10
  •  2
  •   tzot    15 年前

    我会用复数做向量。我的例子是python,它有内置的复数:

    import cmath # complex math
    
    def average_angle(list_of_angles):
    
        # make a new list of vectors
        vectors= [cmath.rect(1, angle) # length 1 for each vector
            for angle in list_of_angles]
    
        vector_sum= sum(vectors)
    
        # no need to average, we don't care for the modulus
        return cmath.phase(vector_sum)
    

    注意,python没有 需要 为了构建一个临时的向量列表,上述所有操作都可以在一个步骤中完成;我只是选择了这种方法来近似适用于其他语言的伪代码。

        11
  •  2
  •   notadam    10 年前

    这里有一个完整的C++解决方案:

    #include <vector>
    #include <cmath>
    
    double dAngleAvg(const vector<double>& angles) {
        auto avgSin = double{ 0.0 };
        auto avgCos = double{ 0.0 };
        static const auto conv      = double{ 0.01745329251994 }; // PI / 180
        static const auto i_conv    = double{ 57.2957795130823 }; // 180 / PI
        for (const auto& theta : angles) {
            avgSin += sin(theta*conv);
            avgCos += cos(theta*conv);
        }
        avgSin /= (double)angles.size();
        avgCos /= (double)angles.size();
        auto ret = double{ 90.0 - atan2(avgCos, avgSin) * i_conv };
        if (ret<0.0) ret += 360.0;
        return fmod(ret, 360.0);
    }
    

    它采用双向量形式的角度,并将平均值简单地作为双向量返回。角度必须以度为单位,当然平均值也以度为单位。

        12
  •  2
  •   Brad Saund    8 年前

    在python中,角度介于[-180,180)

    def add_angles(a, b):
      return (a + b + 180) % 360 - 180
    
    def average_angles(a, b):
      return add_angles(a, add_angles(-a, b)/2)
    

    细节:

    平均 两个角度 有两个平均值相距180度,但我们可能需要更接近的平均值。

    从视觉上看,蓝色的平均值( 和绿色( )生成teal点:

    Original

    角度“环绕”(例如355+10=5),但标准算术将忽略此分支点。 但是,如果角度 与分支点相反,那么( + G )/2给出了最接近的平均值:蓝绿色点。

    对于任意两个角度,我们可以旋转问题,使其中一个角度与分支点相反,执行标准平均,然后再旋转回来。

    rotated returned

        13
  •  1
  •   John with waffle    17 年前

    这里有一个想法:通过反复计算最接近的角度的平均值来建立平均值,保持一个权重。

    另一个想法是:找出给定角度之间的最大间隙。找到平分它的点,然后选取圆上的另一个点作为参考零点,以计算平均值。

        14
  •  1
  •   Zach Scrivena    17 年前

    让我们用圆周长上的点来表示这些角。

    我们能假设所有这些点都落在圆的同一半上吗?(否则,没有明显的方法来定义“平均角度”。想想直径上的两点,例如0度和180度——平均90度还是270度?当我们有3个或更多均匀分布的点时会发生什么?)

    基于这个假设,我们选取半圆上的一个任意点作为“原点”,并测量相对于这个原点的给定角度集(称之为“相对角度”)。请注意,相对角的绝对值严格小于180度。最后,取这些相对角的平均值得到所需的平均角(当然是相对于我们的原点)。

        15
  •  1
  •   cffk    14 年前

    没有单一的“正确答案”。我建议你读这本书, K.V.Mardia和P.E.Jupp,“方向统计”,(Wiley,1999年), 进行彻底的分析。

        16
  •  1
  •   bgp2000    9 年前

    这里是一个使用移动平均值并注意使值正常化的完全算术解。如果所有角度都在圆的一侧(彼此之间的角度在180°以内),则速度很快,并提供正确的答案。

    它在数学上相当于添加偏移量,将值移动到范围(0,180),计算平均值,然后减去偏移量。

    这些注释描述了特定值在任何给定时间都可以具有的范围

    // angles have to be in the range [0, 360) and within 180° of each other.
    // n >= 1
    // returns the circular average of the angles int the range [0, 360).
    double meanAngle(double* angles, int n)
    {
        double average = angles[0];
        for (int i = 1; i<n; i++)
        {
            // average: (0, 360)
            double diff = angles[i]-average;
            // diff: (-540, 540)
    
            if (diff < -180)
                diff += 360;
            else if (diff >= 180)
                diff -= 360;
            // diff: (-180, 180)
    
            average += diff/(i+1);
            // average: (-180, 540)
    
            if (average < 0)
                average += 360;
            else if (average >= 360)
                average -= 360;
            // average: (0, 360)
        }
        return average;
    }
    
        17
  •  1
  •   neilireson    7 年前

    好吧,我参加聚会迟到了,但我想我会加上我的2美分,因为我真的找不到任何确切的答案。最后,我实现了以下的MITSUTA方法的Java版本,我希望提供一种简单而健壮的解决方案。特别是当标准偏差同时提供测量分散度时,如果sd==90,则表示输入角度导致的平均值不明确。

    编辑: 事实上,我意识到,考虑到其他答案中的所有对话和三角法,我最初的实现可以进一步简化,事实上令人担忧的简单。

    /**
     * The Mitsuta method
     *
     * @param angles Angles from 0 - 360
     * @return double array containing
     * 0 - mean
     * 1 - sd: a measure of angular dispersion, in the range [0..360], similar to standard deviation.
     * Note if sd == 90 then the mean can also be its inverse, i.e. 360 == 0, 300 == 60.
     */
    public static double[] getAngleStatsMitsuta(double... angles) {
        double sum = 0;
        double sumsq = 0;
        for (double angle : angles) {
            if (angle >= 180) {
                angle -= 360;
            }
            sum += angle;
            sumsq += angle * angle;
        }
    
        double mean = sum / angles.length;
        return new double[]{mean <= 0 ? 360 + mean: mean, Math.sqrt(sumsq / angles.length - (mean * mean))};
    }
    

    …对于所有你(Java)的怪胎来说,你可以用上面的方法得到一个线的平均角度。

    Arrays.stream(angles).map(angle -> angle<180 ? angle: (angle-360)).sum() / angles.length;
    
        18
  •  0
  •   anubhava    17 年前

    阿尼塔克有正确的解决方案。Nick Fortescue的解决方案在功能上是相同的。

    在哪的特殊情况下

    (sum(x_分量)=0.0&sum(y_分量)=0.0)//例如,2个角度为10。和190。度数ea.

    使用0.0度作为总和

    从计算上讲,您必须测试自atan2(0)以来的这种情况。,0.)未定义,将生成错误。

        19
  •  0
  •   ackb    17 年前

    平均角度phi_avg应具有总和_i_phi_avg-phi_2变为最小的性质,其中差异必须在[-pi,pi)中(因为从另一个方向看可能更短!).通过将所有输入值标准化为[0,2pi),保持运行平均功率因数运行,并选择标准化功率因数运行至[-pi,pi),很容易实现这一点。 (通过加减2pi)。上面的大多数建议都是做其他事情的 具有最小的属性,即平均值 某物 ,但不是角度。

        20
  •  0
  •   ronalchn Damien    13 年前

    英语:

    1. 第二个数据集,所有角度移动180。
    2. 取两个数据集的方差。
    3. 取方差最小的数据集的平均值。
    4. 如果该平均值来自移位集,则将答案再次移动180。

    在蟒蛇中:

    一个numpy nx1角度数组

    if np.var(A) < np.var((A-180)%360):
        average = np.average(A)
    
    else:
        average = (np.average((A-180)%360)+180)%360
    
        21
  •  0
  •   waterworld    13 年前

    (只想分享我从估计理论或统计推断的观点)

    Nimble的试验是得到一组角度的mmse估计值,但它是找到“平均”方向的选择之一;你也可以找到mmae估计值,或者其他一些估计值是“平均”方向,这取决于你的度量量化方向误差;或者更普遍地说,在估计理论中,成本函数的定义。

    ^ mmse/mmae对应最小均方/绝对误差。

    Ackb说,“平均角度phi_avg应该具有总和phi_avg-phi_2变为最小的性质……它们是平均值,但不是角度。”

    ——你用均方意义来量化误差,这是最常见的方法之一,但不是唯一的方法。这里大多数人喜欢的答案(即单位向量和得到结果的角度)实际上是合理的解决方案之一。如果向量的方向被建模为von-mises分布,那么它就是(可以证明)作为我们想要的“平均”方向的ml估计量。这种分布并不复杂,只是一个周期性抽样分布,从二维高斯。见Eqn。(2.179)在主教的书《模式识别和机器学习》中。再说一次,这绝不是唯一代表“平均”方向的最好方法,然而,这是一个既有良好的理论依据又有简单的实施的非常合理的方法。

    Nimble说,“Ackb是对的,这些基于向量的解不能被视为角度的真实平均值,它们只是单位向量对应值的平均值。”

    ——这不是真的。“单位矢量对应物”显示矢量方向的信息。角度是一个不考虑向量长度的量,单位向量是一个附加信息,即长度为1。你可以定义你的“单位”向量为长度2,这并不重要。

        22
  •  0
  •   konsnos    11 年前

    我在@david_hanak的回答帮助下解决了这个问题。 正如他所说:

    当保持在同一个半圆内时,在其他两个半圆之间指向“中间”的角度,例如,对于355和5,这将是0,而不是180。为此,您需要检查两个角度之间的差异是否大于180。如果是这样,在使用上述公式之前,将较小的角度增加360。

    所以我做的是计算所有角度的平均值。然后所有小于这个的角,增加360。然后将它们全部相加,再除以它们的长度,重新计算平均值。

            float angleY = 0f;
            int count = eulerAngles.Count;
    
            for (byte i = 0; i < count; i++)
                angleY += eulerAngles[i].y;
    
            float averageAngle = angleY / count;
    
            angleY = 0f;
            for (byte i = 0; i < count; i++)
            {
                float angle = eulerAngles[i].y;
                if (angle < averageAngle)
                    angle += 360f;
                angleY += angle;
            }
    
            angleY = angleY / count;
    

    工作得很好。

        23
  •  0
  •   E.Rooijen    10 年前

    python函数:

    from math import sin,cos,atan2,pi
    import numpy as np
    def meanangle(angles,weights=0,setting='degrees'):
        '''computes the mean angle'''
        if weights==0:
             weights=np.ones(len(angles))
        sumsin=0
        sumcos=0
        if setting=='degrees':
            angles=np.array(angles)*pi/180
        for i in range(len(angles)):
            sumsin+=weights[i]/sum(weights)*sin(angles[i])
            sumcos+=weights[i]/sum(weights)*cos(angles[i])
        average=atan2(sumsin,sumcos)
        if setting=='degrees':
            average=average*180/pi
        return average
    
        24
  •  0
  •   Martin006    10 年前

    您可以在matlab中使用此函数:

    function retVal=DegreeAngleMean(x) 
    
    len=length(x);
    
    sum1=0; 
    sum2=0; 
    
    count1=0;
    count2=0; 
    
    for i=1:len 
       if x(i)<180 
           sum1=sum1+x(i); 
           count1=count1+1; 
       else 
           sum2=sum2+x(i); 
           count2=count2+1; 
       end 
    end 
    
    if (count1>0) 
         k1=sum1/count1; 
    end 
    
    if (count2>0) 
         k2=sum2/count2; 
    end 
    
    if count1>0 && count2>0 
       if(k2-k1 >= 180) 
           retVal = ((sum1+sum2)-count2*360)/len; 
       else 
           retVal = (sum1+sum2)/len; 
       end 
    elseif count1>0 
        retVal = k1; 
    else 
        retVal = k2; 
    end 
    
        25
  •  0
  •   Ginés Hidalgo    9 年前

    对于任何编程语言,您都可以在以下链接中看到解决方案和一些解释: https://rosettacode.org/wiki/Averages/Mean_angle

    例如, C++解决方案 :

    #include<math.h>
    #include<stdio.h>
    
    double
    meanAngle (double *angles, int size)
    {
      double y_part = 0, x_part = 0;
      int i;
    
      for (i = 0; i < size; i++)
        {
          x_part += cos (angles[i] * M_PI / 180);
          y_part += sin (angles[i] * M_PI / 180);
        }
    
      return atan2 (y_part / size, x_part / size) * 180 / M_PI;
    }
    
    int
    main ()
    {
      double angleSet1[] = { 350, 10 };
      double angleSet2[] = { 90, 180, 270, 360};
      double angleSet3[] = { 10, 20, 30};
    
      printf ("\nMean Angle for 1st set : %lf degrees", meanAngle (angleSet1, 2));
      printf ("\nMean Angle for 2nd set : %lf degrees", meanAngle (angleSet2, 4));
      printf ("\nMean Angle for 3rd set : %lf degrees\n", meanAngle (angleSet3, 3));
      return 0;
    }
    

    输出:

    Mean Angle for 1st set : -0.000000 degrees
    Mean Angle for 2nd set : -90.000000 degrees
    Mean Angle for 3rd set : 20.000000 degrees
    

    MATLAB解决方案 :

    function u = mean_angle(phi)
        u = angle(mean(exp(i*pi*phi/180)))*180/pi;
    end
    
     mean_angle([350, 10])
    ans = -2.7452e-14
     mean_angle([90, 180, 270, 360])
    ans = -90
     mean_angle([10, 20, 30])
    ans =  20.000
    
        26
  •  0
  •   John Thoits    9 年前

    虽然Starblue的答案给出了平均单位向量的角度,但如果您接受在0到2*pi(或0_°到360_°范围内可能有多个答案,则可以将算术平均值的概念扩展到角度。例如,0度和180度的平均值可以是90度或270度。

    算术平均值具有单值和输入值的平方距离之和最小的特性。两个单位向量沿单位圆的距离可以很容易地计算为其点积的反余弦。如果我们通过最小化向量和每个输入单位向量的点积的平方反余弦之和来选择单位向量,那么我们就有一个等价的平均值。同样,记住在特殊情况下可能有两个或两个以上的最小值。

    这个概念可以扩展到任何尺寸,因为沿着单位球面的距离可以用与沿着单位圆的距离(两个单位向量的点积的反余弦)完全相同的方式计算。

    对于圆,我们可以用多种方法求出这个平均值,但我提出了以下O(n^2)算法(角度以弧度表示,我避免计算单位向量):

    var bestAverage = -1
    double minimumSquareDistance
    for each a1 in input
        var sumA = 0;
        for each a2 in input
            var a = (a2 - a1) mod (2*pi) + a1
            sumA += a
        end for
        var averageHere = sumA / input.count
        var sumSqDistHere = 0
        for each a2 in input
            var dist = (a2 - averageHere + pi) mod (2*pi) - pi // keep within range of -pi to pi
            sumSqDistHere += dist * dist
        end for
        if (bestAverage < 0 OR sumSqDistHere < minimumSquareDistance) // for exceptional cases, sumSqDistHere may be equal to minimumSquareDistance at least once. In these cases we will only find one of the averages
            minimumSquareDistance = sumSqDistHere
            bestAverage = averageHere
        end if
    end for
    return bestAverage
    

    如果所有的角度都在彼此180°之内,那么我们可以使用一个更简单的O(n)+O(排序)算法(同样使用弧度并避免使用单位向量):

    sort(input)
    var largestGapEnd = input[0]
    var largestGapSize = (input[0] - input[input.count-1]) mod (2*pi)
    for (int i = 1; i < input.count; ++i)
        var gapSize = (input[i] - input[i - 1]) mod (2*pi)
        if (largestGapEnd < 0 OR gapSize > largestGapSize)
            largestGapSize = gapSize
            largestGapEnd = input[i]
        end if
    end for
    double sum = 0
    for each angle in input
        var a2 = (angle - largestGapEnd) mod (2*pi) + largestGapEnd
        sum += a2
    end for
    return sum / input.count
    

    要使用度数,只需将pi替换为180。如果您计划使用更多维度,那么您很可能需要使用迭代方法来求解平均值。

        27
  •  0
  •   Community Mohan Dere    9 年前

    基于 Alnitak's answer 我编写了一个计算多角度平均值的Java方法:

    如果你的角度是弧度:

    public static double averageAngleRadians(double... angles) {
        double x = 0;
        double y = 0;
        for (double a : angles) {
            x += Math.cos(a);
            y += Math.sin(a);
        }
    
        return Math.atan2(y, x);
    }
    

    如果角度为度:

    public static double averageAngleDegrees(double... angles) {
        double x = 0;
        double y = 0;
        for (double a : angles) {
            x += Math.cos(Math.toRadians(a));
            y += Math.sin(Math.toRadians(a));
        }
    
        return Math.toDegrees(Math.atan2(y, x));
    }
    
        28
  •  0
  •   DynamicChart    8 年前

    问题非常简单。 1。确保所有角度都在-180和180度之间。 2。a将所有非负角度相加,取其平均值,并计算出多少。 2。b.加上所有的负角度,取它们的平均值并计算出多少。 三。取正平均值减去负平均值的差 如果差异大于180,则将差异更改为360减去差异。否则只需改变差异的符号。请注意,差异总是非负的。 平均角等于位置平均值加上差乘以“重量”,负计数除以负计数和正计数之和。

        29
  •  0
  •   Robert Sutton    8 年前

    这里有一些用于平均角度的Java代码,我认为它是相当健壮的。

    public static double getAverageAngle(List<Double> angles)
    {
        // r = right (0 to 180 degrees)
    
        // l = left (180 to 360 degrees)
    
        double rTotal = 0;
        double lTotal = 0;
        double rCtr = 0;
        double lCtr = 0;
    
        for (Double angle : angles)
        {
            double norm = normalize(angle);
            if (norm >= 180)
            {
                lTotal += norm;
                lCtr++;
            } else
            {
                rTotal += norm;
                rCtr++;
            }
        }
    
        double rAvg = rTotal / Math.max(rCtr, 1.0);
        double lAvg = lTotal / Math.max(lCtr, 1.0);
    
        if (rAvg > lAvg + 180)
        {
            lAvg += 360;
        }
        if (lAvg > rAvg + 180)
        {
            rAvg += 360;
        }
    
        double rPortion = rAvg * (rCtr / (rCtr + lCtr));
        double lPortion = lAvg * (lCtr / (lCtr + rCtr));
        return normalize(rPortion + lPortion);
    }
    
    public static double normalize(double angle)
    {
        double result = angle;
        if (angle >= 360)
        {
            result = angle % 360;
        }
        if (angle < 0)
        {
            result = 360 + (angle % 360);
        }
        return result;
    }
    
        30
  •  -3
  •   Rich Bradshaw    15 年前

    我有一个不同于@starblue的方法,它对上面给出的一些角度给出了“正确”的答案。例如:

    • 角度平均值([350,10])=0
    • 角度平均值([-90,90,40])=13.333
    • 角度平均值([350,2])=356

    它对连续角度之间的差异求和。 代码(在matlab中):

    function [avg] = angle_avg(angles)
    last = angles(1);
    sum = angles(1);
    for i=2:length(angles)
        diff = mod(angles(i)-angles(i-1)+ 180,360)-180
        last = last + diff;
        sum = sum + last;
    end
    avg = mod(sum/length(angles), 360);
    end