2017-09-12 60 views
0

在數學我可以轉換多變量矩在累積量和回用MomentConvert:多變量累積量和矩在Python

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

爲一個可以wolframcloud嘗試。

我想這樣做完全一樣的蟒蛇。是否有任何的Python庫這個能力?

+0

你有什麼到目前爲止已經試過?附上一些代碼[最小,完整和可驗證示例](https://stackoverflow.com/help/mcve),並使用[提供的降價選項](https://stackoverflow.com/editing-幫幫我)。 –

+0

那麼我找到了https://programtalk.com/python-examples/statsmodels.stats.moment_helpers./ 但是這隻適用於一個變量。我需要一個庫來做這件事。 – manifold

回答

0

至少一個方向我現在自己編程:

# 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。變量

指數可以通過重複變量來表示。