2015-08-21 9 views
3

У меня есть экспериментальные данные форм (X,Y) и теоретическая модель формы (x(t;*params),y(t;*params)) где t является физическим (но ненаблюдаемым) переменным, и *params являются параметрами, которые я хочу, чтобы определить. t - это непрерывная переменная, и между моделью x и t существует соотношение 1: 1 и между y и t.Подгонка параметрических кривые в Python

В идеальном мире я бы знал значение T (реальное значение параметра) и смог бы выполнить предельно базовые наименьшие квадраты, чтобы найти значения *params. (Обратите внимание, что я не пытается «соединить» значение x и y в моем участке, как и в 31243002 или 31464345.) Я не могу гарантировать, что в моих реальных данных, скрытое значение T монотонно, так как мои данные собранных через несколько циклов.

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

  1. Выберите некоторое значение *params и применить его к модели
  2. Возьмет массив t значений и поместить его в модель, чтобы создать массив model(*params) = (x(*params),y(*params))
  3. интерполировать X (значения данных) в model получить Y_predicted
  4. Выполнить наименьших квадратов (или другие) сравнение между Y и Y_predicted
  5. ли это получить для нового набора *params
  6. В конце концов, выбрать лучшие значения для *params

Есть несколько очевидных проблем с этим подходом.

1) Я недостаточно опытен с кодированием, чтобы разработать очень хорошее «сделать это снова», кроме «попробовать все в пространстве решений», возможно, «попробуй все в грубой сетке», а затем «попробуй все» снова в немного более тонкой сетке в горячих точках грубой сетки ». Я пробовал использовать методы MCMC, но я никогда не находил никаких оптимальных значений, в основном из-за проблемы 2

2) Шаги 2-4 совершенно неэффективны сами по себе.

Я пробовал что-то вроде (наподобие псевдокода, выполняются фактические функции). Существует много незначительных мелочей, которые могут быть связаны с использованием широковещательной передачи на A, B, но они менее значительны, чем проблема, требующая интерполяции для каждого отдельного шага.

Люди, которых я знаю, рекомендовали использовать какой-то алгоритм Expectation Maximization, но я не знаю достаточно об этом, чтобы закодировать один из них с нуля. Я действительно надеюсь, что есть какой-то потрясающий scipy (или иначе открытый исходный) алгоритм, который я не смог найти, который покрывает всю мою проблему, но на данный момент я не надеюсь.

import numpy as np 
import scipy as sci 
from scipy import interpolate 

X_data 
Y_data 

def x(t,A,B): 
    return A**t + B**t 
def y(t,A,B): 
    return A*t + B 

def interp(A,B): 
    ts = np.arange(-10,10,0.1) 
    xs = x(ts,A,B) 
    ys = y(ts,A,B) 
    f = interpolate.interp1d(xs,ys) 
    return f 

N = 101 
lsqs = np.recarray((N**2),dtype=float) 

count = 0 
for i in range(0,N): 
    A = 0.1*i   #checks A between 0 and 10 
    for j in range(0,N): 
     B = 10 + 0.1*j #checks B between 10 and 20 

     f = interp(A,B) 
     y_fit = f(X_data) 
     squares = np.sum((y_fit - Y_data)**2) 

     lsqs[count] = (A,b,squares) #puts the values in place for comparison later 
     count += 1  #allows us to move to the next cell 

i = np.argmin(lsqs[:,2]) 

A_optimal = lsqs[i][0] 
B_optimal = lsqs[i][1] 
+1

Вы должны проверить 'scipy.interpolate.splprep' – rth

+1

Это, к сожалению, полезно только для подключения точек на мои данные' (X, Y) ', а не для установки, чтобы найти' * params' – Necarion

ответ

0

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

x = a+r cos(t) 
y = b+r sin(t) 

при различных значениях t.

В этом случае то, что я хотел бы сделать это устранить переменную t, чтобы получить соотношение между x и y - в данном случае, (x-a)^2+(y-b)^2 = r^2. Если ваши данные идеально подходят для модели, у вас будет (x-a)^2+(y-b)^2 = r^2 в каждом из ваших точек данных. С какой-то ошибкой, вы можете еще найти (a,b,r) минимизировать Eliminate команды

sum_i ((x_i-a)^2 + (y_i-b)^2 - r^2)^2.

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

PS Вы можете сделать лучше в stats.stackexchange, math.stackexchange или mathoverflow.net. Я знаю, что у последней есть страшная репутация, но мы не кусаемся, действительно!

+0

Спасибо за отзыв , Я пытаюсь сделать именно то, что вы предлагаете, - найти эквиваленты a и b. Однако уравнения не всегда аналитически решаются для t, а когда они это делают, они заканчиваются рядом особенностей. Когда я сохраняю его в терминах 't', у меня есть только особенности в' t = 0'. Однако, когда я устраняю 't', я получаю части уравнений, которые выглядят (некоторые константы выведены для простоты) ' y - x + c '= y/x + (y/(xy))^(1/3) ' Это означает, что относительно небольшие ошибки в' x' могут подавить процедуру подгонки. – Necarion

 Смежные вопросы

  • Нет связанных вопросов^_^