Как вы вычисляете среднее значение набора циклических данных?
Я хочу рассчитать среднее значение набора круговых данных. Например, у меня может быть несколько образцов из чтения компаса. Проблема, конечно, в том, как справиться с оберткой. Тот же алгоритм может быть полезен для циферблата.
фактический вопрос сложнее - что означает статистика на сфере или в алгебраическом пространстве, которое "обертывается", например, аддитивная группа mod n. ответ может быть не единственным, например, среднее значение 359 градусов и 1 степени может быть 0 градусов или 180, но статистически 0 выглядит лучше.
Это реальная проблема программирования для меня, и я пытаюсь сделать так, чтобы это не выглядело как просто математическая проблема.
29 ответов:
этот вопрос подробно рассматривается в книге: "Статистика по сферам", Джеффри С. Уотсон, Университет Арканзаса лекция Примечания в математических науках, 1983 John Wiley & Sons, Inc. как уже упоминалось в http://catless.ncl.ac.uk/Risks/7.44.html#subj4 Брюс Карш.
хороший способ оценить средний угол, A, из набора угловых измерений a[i] 0
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, а затем идут против часовой стрелки. Математика должна работать одинаково независимо.
ДЛЯ ЧАСТНОГО СЛУЧАЯ ДВУХ УГЛОВ:
ответ ((a + b) mod 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
ackb правильно, что эти векторные решения не могут считаться истинными средними углами, они являются только средним значением единичных векторных аналогов. Однако предложенное ackb решение не кажется математически обоснованным.
ниже приводится решение, которое математически выводится из цели минимизации (угол[i] - средний угол)^2 (где разница корректируется при необходимости), что делает его истинным арифметическим средним углов.
во-первых, мы нужно посмотреть, в каких именно случаях разница между углами отличается от разницы между их нормальным числом аналогов. Рассмотрим углы x и y, если y >= x - 180 и y
вот изображение, демонстрирующее, где возникают проблемы при вычислении разности углов. Если x лежит в серой области, то будет a проблема.
чтобы минимизировать переменную, в зависимости от кривой, мы можем взять производную от того, что мы хотим минимизировать, а затем найти точку поворота (где производная = 0).
здесь мы применим идею минимизации квадрата разности для получения общей формулы среднего арифметического: sum (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] - x c = подмножество a, где правильная (угловая) разность (a[i]-360) - x cn = размер c d = подмножество 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
одного этого недостаточно, чтобы получить минимум, в то время как он работает для нормальных значений, которые имеют неограниченный набор, поэтому результат определенно будет лежать в диапазоне набора и поэтому действителен. Мы нужен минимум в пределах диапазона (определяемого сегментом). Если минимум меньше нижней границы нашего сегмента, то минимум этого сегмента должен быть на нижней границе (поскольку квадратичные кривые имеют только 1 точку поворота), а если минимум больше верхней границы нашего сегмента, то минимум сегмента находится на верхней границе. После того, как у нас есть минимум для каждого сегмента, мы просто находим тот, который имеет наименьшее значение для того, что мы минимизируем (sum ((b[i]-x)^2) + sum (((c[i]-360) - b)^2) + sum (((d[i]+360)-c)^2)).
вот изображение кривой, которое показывает, как она изменяется в точках, где x=(a[i]+180)%360. Речь идет о наборе данных {65,92,230,320,250}.
вот реализация алгоритма в Java, включая некоторые оптимизации, его сложность O (nlogn). Он может быть уменьшен до O (n), если вы замените сортировку на основе сравнения на сортировку без сравнения, такую как radix род.
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, однако хорошо известно, что среднее арифметическое сильно зависит от значений ребер. Вы также должны помнить, что углы не являются векторами, как это может показаться при работе с углами иногда.
этот алгоритм, конечно, также применяется ко всем величинам, которые подчиняются модульной арифметике (с минимальной корректировкой), например, времени суток.
Я также хотел бы подчеркнуть, что, хотя это истинное среднее значение углов, в отличие от векторных решений, это не обязательно означает, что это решение, которое вы должны использовать, среднее значение соответствующих единичных векторов вполне может быть значением, которое вы на самом деле должны использовать.
вы должны определить в среднем более точно. Для конкретного случая двух углов я могу придумать два разных сценария:
- " истинное " среднее, т. е. (a + b) / 2% 360.
- угол, который указывает "между" двумя другими, оставаясь в том же полукруге, например, для 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 on Matrix Analysis and Applications, Volume 24, Issue 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 сырых показаний подшипников, чтобы сгладить вариации.
- проверьте, является ли первый подшипник диапазоном 270-360 или 0-90 градусов (северные два квадранта)
- если это так, поверните это и все последующие показания на 180 градусов, сохраняя все значения в диапазоне 0
- после того, как 10 показаний были приняты вычислить численное среднее, предполагая, что не было никакого обтекания
- если вращение на 180 градусов было в действительности, то поверните рассчитанное среднее значение на 180 градусов, чтобы вернуться к "истинному" подшипнику.
Это не идеально; он может сломаться. В этом случае мне это сошло с рук, потому что устройство вращается очень медленно. Я поставлю его там на случай, если кто-то еще будет работать под подобное ограничение.
вот полное решение: (входной сигнал массив подшипника в градусах (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; }
Я бы пошел векторным путем, используя комплексные числа. Мой пример на 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); }
Он принимает углы в виде вектора двойников и возвращает среднее значение просто как двойник. Угол должен быть в градусах, и конечно же средний в степени.
в 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° друг от друга, но мы можем ближе средний.
визуально, среднее значение синего (b) и зеленый (a) дает чирок точка:
углы ' обернуть вокруг '(например, 355 + 10 = 5), но стандартная арифметика будет игнорировать это точка ветвления. Однако если угол b находится напротив точки ветвления, то (b + g) / 2 дает самое близкое среднее: чирок точка.
для любых двух углов мы можем повернуть задачу так, чтобы один из углов был противоположен точке ветвления, выполнить стандартное усреднение, а затем повернуть назад.
вот идея: построить среднее итеративно, всегда вычисляя среднее из углов, которые находятся ближе всего друг к другу, сохраняя вес.
еще одна идея: найти наибольший зазор между заданными углами. Найдите точку, которая делит его пополам, а затем выберите противоположную точку на окружности в качестве опорного нуля для вычисления среднего значения.
обозначим эти углы с точками на окружности.
можно ли предположить, что все эти точки попадают на одну и ту же половину окружности? (В противном случае нет очевидного способа определить "средний угол". Подумайте о двух точках на диаметре, например 0 град и 180 град - - - это средний 90 град или 270 град? Что происходит, когда у нас есть 3 или более равномерного распределения точек?)
с этим предположением, мы выбираем произвольную точку на полуокружности как "начало координат "и измеряют заданный набор углов относительно этого начала координат (назовем это"относительным углом"). Обратите внимание, что относительный угол имеет абсолютное значение строго меньше 180 градусов. Наконец, возьмите среднее значение этих относительных углов, чтобы получить желаемый средний угол (относительно нашего происхождения, конечно).
нет ни одного"правильного ответа". Рекомендую почитать книгу, К. В. Мардия и П. Е. Юпп, "направленная статистика", (Wiley, 1999), для тщательного анализа.
Альнитак имеет правильное решение. Решение Ника Фортескью функционально то же самое.
для частного случая где
(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| to [- Pi,Pi) (путем сложения или вычитания 2Pi). Большинство предложений выше сделать что-то еще, что делает не имейте это минимальное свойство, т. е. в среднем что-то, но не углы.
на английском языке:
- сделайте второй набор данных со всеми углами, сдвинутыми на 180.
- возьмите дисперсию обоих наборов данных.
- возьмите среднее значение набора данных с наименьшей дисперсией.
- если это среднее значение из сдвинутого набора, то снова сдвиньте ответ на 180.
в python:
массив углов #numpy NX1
if np.var(A) < np.var((A-180)%360): average = np.average(A) else: average = (np.average((A-180)%360)+180)%360
(просто хочу поделиться своей точкой зрения из теории оценки или статистического вывода)
испытание Nimble заключается в том, чтобы получить оценку MMSE^ набора углов, но это один из вариантов, чтобы найти "усредненное" направление; можно также найти оценку MMAE^ или какую-либо другую оценку, чтобы быть "усредненным" направлением, и это зависит от вашей метрической количественной ошибки направления; или, более широко, в теории оценки, определение функции стоимости.
^ MMSE / MMAE соответствует к минимальной среднеквадратичной / абсолютной ошибке.
ackb сказал: "средний угол phi_avg должен иметь свойство, что sum_i|phi_avg-phi_i / ^2 становится минимальным...они усредняют что-то, но не углы"
---- вы количественно оцениваете ошибки в среднеквадратичном смысле, и это один из наиболее распространенных способов, однако, не единственный способ. Ответ, который предпочитают большинство людей здесь (т. е. сумма единичных векторов и получить угол результата), на самом деле является одним из разумных решений. Это (может быть доказано) оценка ML, которая служит "усредненным" направлением, которое мы хотим, если направления векторов моделируются как распределение фон Мизеса. Это распределение не фантазии, а просто периодически выборочные распределения из 2D Guassian. См. Eqn. (2.179) в книге Бишопа "распознавание образов и машинное обучение". Опять же, ни в коем случае это не единственный лучший способ представить "среднее" направление, однако вполне разумный, который имеет как хорошее теоретическое обоснование, так и простое реализация.
Nimble сказал: "ackb прав, что эти векторные решения не могут считаться истинными средними углами, они являются только средним значением единичных векторных аналогов"
----Это неправда. "Единичные векторные аналоги" раскрывают информацию о направлении вектора. Угол-это величина без учета длины вектора, а единичный вектор-это что-то с дополнительной информацией о том, что длина равна 1. Вы можете определить свой "блок" вектор должен быть длиной 2, это действительно не имеет значения.
Я решил проблему с помощью ответа от @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*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) (снова используя радианы и избегая использования единицы измерения векторы):
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. Если вы планируете использовать больше измерений, то вам, скорее всего, придется использовать итерационный метод для решения для среднего.
вот полностью арифметическое решение, используя скользящие средние и заботясь о нормализации значений. Это быстро и дает правильные ответы, если все углы на одной стороне круга (в пределах 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; }
на основе ответ Алнитака, Я написал метод 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)); }
задача предельно проста. 1. Убедитесь, что все углы между -180 и 180 градусов. 2. A добавьте все неотрицательные углы, возьмите их среднее значение и подсчитайте, сколько 2. B. сложите все отрицательные углы, возьмите их среднее значение и посчитайте, сколько их. 3. Возьмите разницу pos_average минус neg_average Если разница больше 180, то измените разницу на 360 минус разница. В противном случае просто измените знак различия. Обратите внимание, что разница всегда неотрицательна. средний угол равняется 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