2017-09-12 29 views
0

に私はキュムラントで多変量瞬間を変換し、バックMomentConvert使用することができます。一つはwolframcloudに試すことができるよう多変数キュムラントとモーメントMathematicaでのpython

MomentConvert[Cumulant[{2, 2,1}], "Moment"] // TraditionalForm 

を。

私はPythonで全く同じことをしたいと思います。これに対応するPythonのライブラリはありますか?

+0

これまでに何を試しましたか? [提供されているマークダウンオプション](https://stackoverflow.com/editing-library/)を使用して、[最小、完全、および検証可能な例](https://stackoverflow.com/help/mcve)助けて)。 –

+0

さて、私はhttps://programtalk.com/python-examples/statsmodels.stats.moment_helpers./ を見つけましたが、これは1つの変数に対してのみです。私はこれを行う図書館が必要です。 – manifold

答えて

0

私は今、自分でプログラムされた少なくとも1つの方向:

# from http://code.activestate.com/recipes/577211-generate-the-partitions-of-a-set-by-index/ 

from collections import defaultdict 

class Partition: 

    def __init__(self, S):   
     self.data = list(S)   
     self.m = len(S) 
     self.table = self.rgf_table() 

    def __getitem__(self, i): 
     #generates set partitions by index 
     if i > len(self) - 1: 
      raise IndexError 
     L = self.unrank_rgf(i) 
     result = self.as_set_partition(L) 
     return result 

    def __len__(self): 
     return self.table[self.m,0] 

    def as_set_partition(self, L): 
     # Transform a restricted growth function into a partition 
     n = max(L[1:]+[1]) 
     m = self.m 
     data = self.data 
     P = [[] for _ in range(n)] 
     for i in range(m): 
      P[L[i+1]-1].append(data[i]) 
     return P 

    def rgf_table(self): 
     # Compute the table values 
     m = self.m 
     D = defaultdict(lambda:1) 
     for i in range(1,m+1): 
      for j in range(0,m-i+1): 
       D[i,j] = j * D[i-1,j] + D[i-1,j+1] 
     return D 

    def unrank_rgf(self, r): 
     # Unrank a restricted growth function 
     m = self.m 
     L = [1 for _ in range(m+1)] 
     j = 1 
     D = self.table 
     for i in range(2,m+1): 
      v = D[m-i,j] 
      cr = j*v 
      if cr <= r: 
       L[i] = j + 1 
       r -= cr 
       j += 1 
      else: 
       L[i] = r // v + 1 
       r %= v 
     return L 



# S = set(range(4)) 
# P = Partition(S) 
# for x in P: 
#  print (x) 


# using https://en.wikipedia.org/wiki/Cumulant#Joint_cumulants 
import math 

def Cum2Mom(arr, state): 
    def E(op): 
     return qu.expect(op, state) 

    def Arr2str(arr,sep): 
     r = '' 
     for i,x in enumerate(arr): 
      r += str(x) 
      if i<len(arr)-1: 
       r += sep     
     return r 

    if isinstance(arr[0],str): 
     myprod = lambda x: Arr2str(x,'*') 
     mysum = lambda x: Arr2str(x,'+') 
     E=lambda x: 'E('+str(x)+')' 
     myfloat = str 
    else: 
     myfloat = lambda x: x 
     myprod = np.prod 
     mysum = sum 
    S = set(range(len(arr))) 
    P = Partition(S) 

    return mysum([ 
     myprod([myfloat(math.factorial(len(pi)-1) * (-1)**(len(pi)-1)) 
     ,myprod([ 
      E(myprod([ 
       arr[i] 
       for i in B 
      ])) 
      for B in pi])]) 
     for pi in P]) 

print(Cum2Mom(['a','b','c','d'],1) ) 
import qutip as qu 
print(Cum2Mom([qu.qeye(3) for i in range(3)],qu.qeye(3)) ) 

qutipのopjectsで動作するように設計されており、それも正しい分離とprefactorsを検証するための文字列を使用していますいます。

変数の指数は、変数を繰り返すことで表すことができます。

関連する問題