На главную Наши проекты:
Журнал   ·   Discuz!ML   ·   Wiki   ·   DRKB   ·   Помощь проекту
ПРАВИЛА FAQ Помощь Участники Календарь Избранное RSS
msm.ru
! правила раздела Алгоритмы
1. Помните, что название темы должно хоть как-то отражать ее содержимое (не создавайте темы с заголовком ПОМОГИТЕ, HELP и т.д.). Злоупотребление заглавными буквами в заголовках тем ЗАПРЕЩЕНО.
2. При создании темы постарайтесь, как можно более точно описать проблему, а не ограничиваться общими понятиями и определениями.
3. Приводимые фрагменты исходного кода старайтесь выделять тегами code.../code
4. Помните, чем подробнее Вы опишете свою проблему, тем быстрее получите вразумительный совет
5. Запрещено поднимать неактуальные темы (ПРИМЕР: запрещено отвечать на вопрос из серии "срочно надо", заданный в 2003 году)
6. И не забывайте о кнопочках TRANSLIT и РУССКАЯ КЛАВИАТУРА, если не можете писать в русской раскладке :)
Модераторы: Akina, shadeofgray
Страницы: (3) [1] 2 3  все  ( Перейти к последнему сообщению )  
> Как определить параметры синуса?
    Есть результаты физического эксперимента, которые с некоторой ошибкой могут быть представлены в виде функции
    y = A*sin(B*x + x0) + y0
    Всего пар x;y около 1000, причем в эксперимент попадает несколько(8-10) периодов функции.
    Иными словами точки насажены довольно плотно.

    Вопрос - каким образом можно определить A, B, x0, y0 так, что бы точность определения параметров росла с ростом числа известных точек, ну и желательно что бы минимизировалась ошибка(но уже необязательно)?
      Это называется "нелинейный МНК". Используй Левенберга-Марквардта для минимизации суммарной ошибки, т.е. функции вида F=(y(x0)-y0)^2 + ... + (y(xn)-yn)^2

      Реализация метода - либо пакет levmar (c/c++), либо из ALGLIB - http://alglib.sources.ru/optimization/levenbergmarquardt.php

      Добавлено
      Добавление 1: да, минимизация ведется по A,B, xo, yo.

      Добавление 2: в такой постановке задачи у тебя минимизируется ошибка аппроксимации (т.е. ошибка вычисления функции y), а не ошибка нахождения неизвестных параметров. Однако на практике можно исходить из того, что малая ошибка аппроксимации обозначает малую ошибку в параметрах функции.
        (спасибо, пока думаю)

        Добавлено
        Тут нужна некоторая помощь.

        У нас выборка x(i), y(i)
        Как я понимаю, в качестве f(i) берется (Y - A*sin(B*x(i) + X0) + Y0)^2
        Где Y,A,X0,Y0 - переменные от которых зависит функция.
        Я взял реализацию на C# и применял ее вот так:
        ExpandedWrap disabled
                  Dim s As New alglib.minlm.lmstate
                  Dim sol() As Double = {0.25, 1591.549431, 0, -150} ' {0.1, 1000, 0, -120} '  'Примерно подобранные начальные значения
                  'ReDim sol(results.Count)
                  alglib.minlm.minlmfj(4, results.Count, sol, 0.0000001, 0.0000001, 10000, s)' Тут числа взяты наугад
           
                  While alglib.minlm.minlmiteration(s)
                      Dim A As Double = s.x(0)
                      Dim B As Double = s.x(1)
                      Dim X0 As Double = s.x(2)
                      Dim Y0 As Double = s.x(3)
                      If s.needf Then
                          Dim _F As Double = 0
                          For Each R As PointF In results
                              _F += (R.Y - (A * Math.Sin(B * R.X + X0) + Y0)) ^ 2
                          Next
                          s.f = _F
                          Debug.Print(_F.ToString) 'ТУТ вывожу значения вычисляемой функции
                      End If
                      If s.needfij Then
                          For i As Integer = 0 To results.Count - 1
                              Dim r As PointF = results(i)
                              s.fi(i) = r.Y - (A * Math.Sin(B * r.X + X0) + Y0)
                              s.j(i, 0) = Math.Sin(B * r.X + X0) REM dF/dA
                              s.j(i, 1) = A * B * Math.Cos(B * r.X + X0) REM dF/dB
                              s.j(i, 2) = A * Math.Cos(B * r.X + X0) REM dF/dX0
                              s.j(i, 3) = 1 REM dF/dY0
                          Next
                      End If
                      If s.needfg Or s.needfgh Then
                          Throw New Exception()
                      End If
                  End While
                  Dim rep As alglib.minlm.lmreport
                  Dim X(4) As Double
                  alglib.minlm.minlmresults(s, X, rep)
                  Dim z As String = "{" & X(0) & ", " & X(1) & ", " & X(2) & ", " & X(3) & "}"
                  MsgBox(z)

        Использую вот такой код.


        Теперь пробую прогнать ее на одной и той же выборке, при разных значениях стартовых параметров
        При стартовых параметрах {0.25, 1591.549431, 0, -150} (это довольно близко к правде) алгоритм дает ответ
        {0.250000043588598, 1591.54946567717, 2.17883057160585E-08, -150.00006958343}
        Что похоже на правду.
        При этом последовательные значения функции вот такие:
        Цитата
        175285943.807378
        175279093.429274
        174998443.565119
        155484743.497379
        47898495.6142789
        43850845.497696


        Если же взять другие начальные значения, например {0.1, 1000, 0, -120} что уже довольно далеко от истины, получается ответ
        {0.0999999859957626, 1000.00000574627, 5.74626982824354E-09, -120.000062631127}
        А это совсем неверно
        При этом вычисленные значения функции следующие:
        Цитата
        142013745.650736
        142008195.060041
        141780793.005068
        125971771.737748
        38809425.2111686
        35530193.4296744

        Вообще какие-то громадные числа получаются. Причина завершения в обоих случаях "relative step is no more than EpsX"
          Ну, насколько я помню этот метод, он ищет локальный минимум, а если их несколько (а в реальных многопараметрических задачах - обычно очень много) - в какой он попадет, как раз зависит от начальных условий.
            Гм. Ну а как тогда все таки решать исходную задачу? Нужно найти, получается, глобальный минимум, причем в нем f(x) = 0
            И посмотрите плз, правильно ли посчитан якобиан
              А нельзя ли для нахождения параметров синуса применить спектральный анализ?
                Что?)
                  Э-э-э... может я не понял суть задачи :unsure: .. в более понятных мне терминах радиотехнических, она звучит так: необходимо определить амплитуду, частоту, начальную фазу и постоянную составляющую синусоидального сигнала, заданного дискретными значениями. Только значения y должны быть "взяты" через равные промежутки времени, а в Вашей формуле вместо времени фигурирует некий x, и о нем ничего неизвестно. Если растояния между соседними значениями x одинаковы, можно вычислить спектр сигнала.
                  Сообщение отредактировано: Prince -
                    Цитата
                    А нельзя ли для нахождения параметров синуса применить спектральный анализ?
                    Да, совершенно верно, я о том же подумал.
                    Левенберг-Марквардт в данном случае - это чересчур.
                      ANDLL, неправильно считается Якобиан.

                      Цитата ANDLL @

                      s.j(i, 0) = Math.Sin(B * r.X + X0) REM dF/dA
                      s.j(i, 1) = A * B * Math.Cos(B * r.X + X0) REM dF/dB
                      s.j(i, 2) = A * Math.Cos(B * r.X + X0) REM dF/dX0
                      s.j(i, 3) = 1 REM dF/dY0


                      Во-первых, в каждой строчке потерян минус. Во-вторых, в строке "s.j(i,1)=" есть лишний множитель B.
                        Гм.
                        Только наверное
                        s.j(i, 1) = - A * r.X * Math.Cos(B * r.X + X0) REM dF/dB
                          Цитата 0797293 @
                          Цитата
                          А нельзя ли для нахождения параметров синуса применить спектральный анализ?

                          Да, совершенно верно, я о том же подумал.
                          Левенберг-Марквардт в данном случае - это чересчур.

                          Мне кажется, это как раз спектральный анализ - не в тему :) Во-первых, методы, основанные на БПФ, позволяют работать только с точками на равномерной сетке - а здесь сетка неравномерная. Во-вторых, спектральные методы позволяют выбирать множитель "B" только из конечного диапазона частот, а Левенберг-Марквардт позволяет непрерывно выбирать оптимальный параметр. В-третьих, спектральный анализ совершенно не годится для ситуации, когда нельзя однозначно определить значение функции в точке (есть несколько точек с совпадающими абсциссами, но разными значениями функции) - ситуация, часто встречающаяся в практике.


                          P.S. Вообще, наверное надо добавить в алгоритм какой-то код для проверки точности вычисления Якобиана/градиента/Гессиана.

                          Добавлено
                          Цитата ANDLL @
                          Гм.
                          Только наверное
                          s.j(i, 1) = - A * r.X * Math.Cos(B * r.X + X0) REM dF/dB

                          А, ну да, конечно же вы правы :) я неясно выразился.
                            Заменил код.
                            Вот так правильный якобиан?
                            ExpandedWrap disabled
                                                  s.j(i, 0) = -Math.Sin(B * r.X + X0) REM dF/dA
                                                  s.j(i, 1) = -A * r.X * Math.Cos(B * r.X + X0) REM dF/dB
                                                  s.j(i, 2) = -A * Math.Cos(B * r.X + X0) REM dF/dX0
                                                  s.j(i, 3) = -1 REM dF/dY0

                            А то результат все равно довольно неправильный(
                              Цитата
                              Во-первых, методы, основанные на БПФ, позволяют работать только с точками на равномерной сетке

                              Вот в этом и была неясность. Раз сетка неравномерная, со спектром я промахнулся.
                                К слову сказать, точки и правда лежат на равномерной сетке(ну или по крайней мере можно выкинув небольшое число точек ее получить) , чето я сразу не подумал что это важно.
                                0 пользователей читают эту тему (0 гостей и 0 скрытых пользователей)
                                0 пользователей:
                                Страницы: (3) [1] 2 3  все


                                Рейтинг@Mail.ru
                                [ Script execution time: 0,0439 ]   [ 15 queries used ]   [ Generated: 17.06.24, 13:09 GMT ]