Как рассчитать среднее значение для набора циклических данных?

Я хочу рассчитать среднее значение для набора циклических данных. Например, у меня может быть несколько образцов из чтения компаса. Проблема, конечно же, заключается в том, как бороться с циклом. Тот же алгоритм может быть полезен для циферблата.

Фактический вопрос более сложен - что означает статистика на сфере или в алгебраическом пространстве, которое «оборачивается», например, аддитивная группа mod n. Ответ может быть не уникальным, например, среднее значение 359 градусов и 1 градус может быть 0 или 180, но статистически 0 выглядит лучше.

Для меня это настоящая проблема программирования, и я пытаюсь сделать так, чтобы она не выглядела просто математической.

Вычислить единичные векторы по углам и взять угол их среднего.

Этот вопрос подробно рассматривается в книге «Статистика по сферам», Джеффри С. Уотсон, Лекционные заметки по математическим наукам Университета Арканзаса, 1983 г., John Wiley & Sons, Inc., как упомянуто по адресу http: //catless.ncl. ac.uk/Risks/7.44.html#subj4 Брюса Карша.

Хороший способ оценить средний угол A из набора угловых измерений a [i] 0 <= i

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

Метод, данный starblue, эквивалентен в вычислительном отношении, но его причины яснее и, вероятно, программно более эффективны, а также хорошо работают в нулевом случае, так что слава ему.

Теперь предмет более подробно рассматривается в Википедии , а также в других областях, таких как дробные части.

Я вижу проблему - например, если у вас угол 45 'и угол 315', «естественное» среднее значение будет 180 ', но на самом деле вы хотите получить значение 0'.

Я думаю, что Starblue на что-то. Просто вычислите (x, y) декартовы координаты для каждого угла и сложите полученные векторы вместе. Угловое смещение конечного вектора должно быть вашим требуемым результатом.

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

Сейчас я игнорирую, что направление по компасу начинается с севера и идет по часовой стрелке, тогда как «нормальные» декартовы координаты начинаются с нуля вдоль оси X, а затем идут против часовой стрелки. Математика должна работать одинаково независимо.

ПО ОСОБЕННОМУ СЛУЧАЮ ДВУХ УГЛОВ:

Ответ ((а + б) мод 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
• 1, 180 -> 90,5
• 20, 350 -> 5
• 350, 20 -> 5 (все последующие примеры тоже меняются местами)
• 10, 20 -> 15
• 350, 2 -> 356
• 359, 0 -> 359,5
• 180, 180 -> 180

Акб прав, что эти векторные решения нельзя считать истинными средними значениями углов, они являются только средним значением единичных векторов. Тем не менее, предложенное решение Акба не кажется математически обоснованным.

Ниже приведено решение, которое математически получено из цели минимизации (angle [i] - avgAngle) ^ 2 (где разница корректируется при необходимости), что делает его истинным средним арифметическим углов.

Во-первых, нам нужно посмотреть, в каких именно случаях разница между углами отличается от разницы между их обычными числами. Рассмотрим углы x и y, если y> = x - 180 и y <= x + 180, то мы можем напрямую использовать разность (xy). В противном случае, если первое условие не выполняется, мы должны использовать (y + 360) в расчете вместо y. Соответственно, если второе условие не выполнено, тогда мы должны использовать (y-360) вместо y. Поскольку в уравнении кривой мы минимизируем только изменения в точках, где эти неравенства изменяются с истинного на ложное или наоборот, мы можем разделить полный диапазон [0,360) на набор сегментов, разделенных этими точками. Затем нам нужно только найти минимум каждого из этих сегментов, а затем минимум каждого сегмента, который является средним.

Вот изображение, демонстрирующее, где возникают проблемы при расчете угловых разностей. Если х лежит в серой области, тогда возникнет проблема.

Чтобы минимизировать переменную, в зависимости от кривой, мы можем взять производную того, что мы хотим минимизировать, и затем мы найдем точку поворота (где производная = 0).

Здесь мы применим идею минимизации квадрата разности для получения общей средней арифметической формулы: сумма (a [i]) / n. Кривая y = sum ((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] -xc = подмножество a, где правильная (угловая) разница (a [i] -360) -x cn = размер cd = подмножество a, где правильная (угловая) разница (a [i] +360) -x dn = размер 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

Одного этого недостаточно для получения минимума, в то время как он работает для нормальных значений, у которых есть неограниченный набор, поэтому результат определенно будет лежать в пределах диапазона набора и, следовательно, действителен. Нам нужен минимум в пределах диапазона (определяемого сегментом). Если минимум меньше нижней границы нашего сегмента, то минимум этого сегмента должен быть на нижней границе (поскольку у квадратичных кривых есть только одна точка поворота), а если минимум больше верхней границы нашего сегмента, то минимум сегмента находится на верхняя граница. После того, как у нас есть минимум для каждого сегмента, мы просто находим тот, который имеет наименьшее значение для того, что мы минимизируем (sum ((b [i] -x) ^ 2) + sum (((c [i] -360 ) -b) ^ 2) + сумма (((d [i] +360) -c) ^ 2)).

Вот изображение кривой, которая показывает, как она изменяется в точках, где x = (a [i] +180)% 360. Набор данных находится под вопросом {65,92,230,320,250}.

Вот реализация алгоритма на 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_;
}

Среднее арифметическое набора углов может не совпадать с вашим интуитивным представлением о том, каким должно быть среднее значение. Например, среднее арифметическое для набора {179,179,0,181,181} равно 216 (и 144). Ответ, о котором вы сразу подумаете, вероятно, равен 180, однако хорошо известно, что среднее арифметическое сильно зависит от значений ребер. Вы также должны помнить, что углы не являются векторами, как это может показаться привлекательным, когда иногда приходится иметь дело с углами.

Этот алгоритм, конечно, также применим ко всем величинам, которые подчиняются модульной арифметике (с минимальной корректировкой), например, ко времени суток.

Я также хотел бы подчеркнуть, что, хотя это истинное среднее углов, в отличие от векторных решений, это не обязательно означает, что это решение, которое вы должны использовать, среднее значение соответствующих единичных векторов вполне может быть значением, которое вы на самом деле следует использовать.

Вы должны определить среднее значение более точно. Для конкретного случая двух углов я могу представить два разных сценария:

1. «Истинное» среднее значение, т.е. (a + b) / 2% 360.
2. Угол, который указывает «между» двумя другими, оставаясь в том же полукруге, например, для 355 и 5, это будет 0, а не 180. Чтобы сделать это, вам нужно проверить, больше ли разница между двумя углами, чем 180 или не. Если это так, увеличьте меньший угол на 360, прежде чем использовать приведенную выше формулу.

Однако я не понимаю, как можно обобщить вторую альтернативу на случай более двух углов.

Как и все средние, ответ зависит от выбора метрики. Для данной метрики M среднее значение некоторых углов a_k в [-pi, pi] для k в [1, N] является тем углом 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)

Два общих выбора метрики - это метрики Фробениуса и Римана. Для метрики Фробениуса существует прямая формула, которая соответствует обычному понятию среднего значения в круговой статистике. См. «Средства и усреднение в группе ротаций», Maher Moakher, SIAM Journal по матричному анализу и приложениям, том 24, выпуск 1, 2002, для получения более подробной информации.
http://link.aip.org/link/?SJMAEL/24/1/1

Вот функция для GNU Octave 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/'.

Я хотел бы поделиться методом, который я использовал с микроконтроллером, который не имел возможностей с плавающей запятой или тригонометрии. Мне все еще нужно было «усреднить» 10 необработанных показаний подшипников, чтобы сгладить вариации.

1. Проверьте, находится ли первый подшипник в диапазоне 270-360 или 0-90 градусов (два северных квадранта)
2. Если это так, поверните это и все последующие показания на 180 градусов, сохраняя все значения в диапазоне 0 <= подшипник <360. В противном случае считайте показания по мере их поступления.
3. После того, как были сделаны 10 показаний, вычислите среднее числовое значение, предполагая, что не было никакого обхода
4. Если вращение на 180 градусов действовало, то поверните рассчитанное среднее на 180 градусов, чтобы вернуться к «истинному» азимуту.

Это не идеально; это может сломаться. В этом случае мне это сошло с рук, потому что устройство вращается очень медленно. Я опубликую это на тот случай, если кто-то окажется, что работает в подобных условиях.

По-английски:

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

Вот полное решение: (входной сигнал представляет собой массив азимутов в градусах (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;
}

В питоне, с углами между [-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 °, но мы можем хотеть более близкое среднее.

Визуально, среднее значение синего ( b ) и зеленого ( a ) значений дает точку размытия:

Углы «окружают» (например, 355 + 10 = 5), но стандартная арифметика игнорирует эту точку ветвления. Однако, если угол b противоположен точке ветвления, то ( b + g ) / 2 дает ближайший средний показатель: чирокую точку.

Для любых двух углов мы можем повернуть задачу так, чтобы один из углов был противоположен точке ветвления, выполнить стандартное усреднение, а затем повернуть назад.

Я бы пошел вектор, используя комплексные числа. Мой пример в 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 не должен создавать временный новый список векторов, все вышеперечисленное можно сделать за один шаг; Я просто выбрал этот способ, чтобы приблизить псевдокод, применимый и к другим языкам.

Вот полное решение 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);
}

Он принимает углы в виде вектора двойных значений и возвращает среднее значение просто как двойной. Углы должны быть в градусах, и, конечно, среднее значение также в градусах.

Основываясь на ответе Алнитак , я написал метод 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));
}

Вот идея: построить среднее значение итеративно, всегда вычисляя среднее из углов, которые находятся ближе всего друг к другу, сохраняя вес.

Другая идея: найти самый большой зазор между заданными углами. Найдите точку, которая делит ее пополам, а затем выберите противоположную точку на круге в качестве контрольного нуля, чтобы вычислить среднее значение.

Представим эти углы точками на окружности круга.

Можно ли предположить, что все эти точки попадают на одну и ту же половину круга? (В противном случае, нет очевидного способа определить «средний угол». Подумайте о двух точках на диаметре, например, 0 и 180 градусов - это в среднем 90 или 270 градусов? Что происходит, когда у нас есть 3 или более равномерно распределить точки?)

С этим предположением мы выбираем произвольную точку на этом полукруге в качестве «начала координат» и измеряем заданный набор углов относительно этого начала координат (назовем это «относительным углом»). Обратите внимание, что относительный угол имеет абсолютное значение строго менее 180 градусов. Наконец, возьмите среднее значение этих относительных углов, чтобы получить желаемый средний угол (относительно нашего происхождения, конечно).

Там нет ни одного "правильного ответа". Я рекомендую прочитать книгу KV Mardia и PE Jupp, "Direction Statistics" (Wiley, 1999), для тщательного анализа.

(Просто хочу поделиться своей точкой зрения из теории оценивания или статистического вывода)

Задача Nimble состоит в том, чтобы получить оценку MMSE ^ для набора углов, но это один из вариантов поиска «усредненного» направления; Можно также найти оценку MMAE ^ или некоторую другую оценку, чтобы быть «усредненным» направлением, и это зависит от вашей количественной погрешности измерения метрики направления; или, в более общем смысле, в теории оценки - определение функции стоимости.

^ MMSE / MMAE соответствует минимальному среднему квадрату / абсолютной ошибке.

Акб сказал: «Средний угол phi_avg должен иметь свойство, которое sum_i | phi_avg-phi_i | ^ 2 становится минимальным ... они что-то усредняют, но не углы»

---- вы определяете количество ошибок в среднеквадратичном смысле, и это один из наиболее распространенных способов, однако не единственный. Ответ, одобренный большинством людей здесь (т. Е. Сумма векторов единиц и получение угла результата), на самом деле является одним из разумных решений. Это (может быть доказано) оценка ML, которая служит «усредненным» направлением, которое мы хотим, если направления векторов моделируются как распределение фон Мизеса. Это распределение не является причудливым и представляет собой просто периодически выбираемое распределение из двумерного гассиана. Смотрите уравнение (2.179) в книге Бишопа «Распознавание образов и машинное обучение». Опять же, ни в коем случае это не единственный лучший способ представить «среднее» направление, однако это вполне разумный вариант, имеющий хорошее теоретическое обоснование и простую реализацию.

Нимбл сказал: «Акб прав, что эти векторные решения нельзя считать истинными средними значениями углов, они являются только средними значениями единичных векторов».

----это неправда. «Единичные векторные аналоги» раскрывают информацию о направлении вектора. Угол - это величина без учета длины вектора, а единичный вектор - это нечто с дополнительной информацией о том, что длина равна 1. Вы можете определить, что ваш «единичный» вектор имеет длину 2, это на самом деле не имеет значения.

Вот полностью арифметическое решение с использованием скользящих средних и заботой о нормализации значений. Это быстро и дает правильные ответы, если все углы находятся на одной стороне круга (в пределах 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;
}

Ну, я очень опоздал на вечеринку, но подумал, что добавлю свои 2 цента, потому что я не смог найти окончательного ответа. В конце я реализовал следующую Java-версию метода Mitsuta, которая, я надеюсь, обеспечивает простое и надежное решение. В частности, стандартное отклонение обеспечивает как дисперсию измерения, так и, если 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;

Альнитак имеет правильное решение. Решение Ника Фортескью функционально одинаково.

Для особого случая, где

(sum (x_component) = 0.0 && sum (y_component) = 0.0) // например, 2 угла 10 и 190. градусов.

используйте 0,0 градусов в качестве суммы

В вычислительном отношении вы должны проверить этот случай, так как atan2 (0, 0) не определено и приведет к ошибке.

Средний угол phi_avg должен обладать свойством того, что sum_i | phi_avg-phi_i | ^ 2 становится минимальным, где разница должна быть в [-Pi, Pi) (потому что это может быть короче, если будет наоборот!). Это легко достигается путем нормализации всех входных значений до [0, 2Pi), сохранения среднего значения phi_run и выбора нормализации | phi_i-phi_run | к [-Pi, Pi) (путем добавления или вычитания 2Pi). Большинство предложений выше делают что-то еще, что не имеет этого минимального свойства, то есть они усредняют что-то , но не углы.

Я решил проблему с помощью ответа от @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;

Работает отлично.

Функция 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

Вы можете использовать эту функцию в 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 

Вы можете увидеть решение и небольшое объяснение в следующей ссылке для ЛЮБОГО языка программирования: 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

В то время как ответ starblue дает угол среднего единичного вектора, можно расширить понятие среднего арифметического на углы, если принять, что может быть более одного ответа в диапазоне от 0 до 2 * пи (или от 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) (снова используя радианы и избегая использования единичных векторов):

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. Если вы планируете использовать больше измерений, вам, скорее всего, придется использовать итерационный метод для вычисления среднего значения.

Проблема чрезвычайно проста. 1. Убедитесь, что все углы находятся в диапазоне от -180 до 180 градусов. 2. a Добавьте все неотрицательные углы, возьмите их среднее и СЧИТАЙТЕ сколько 2. б. Добавьте все отрицательные углы, возьмите их среднее и СЧИТАЙТЕ сколько. 3. Возьмите разницу pos_average минус neg_average. Если разница больше 180, то измените разницу на 360 минус разница. В противном случае просто измените знак различия. Обратите внимание, что разница всегда неотрицательна. Average_Angle равен pos_average плюс разница, умноженная на «вес», отрицательное число, деленное на сумму отрицательного и положительного числа

Вот немного кода 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;
}

У меня есть метод, отличный от @Starblue, который дает «правильные» ответы на некоторые из приведенных выше углов. Например:

• angle_avg ([350,10]) = 0
• angle_avg ([- 90,90,40]) = 13.333
• angle_avg ([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