"""
Function of this code: one method of feature extraciton from the numerical values.
Contact for help: gch_liu@163.com 
"""
import numpy as np

def GetMomentVectorPS(seq):
    n=len(seq);
    uA=np.zeros(n);
    uC=np.zeros(n);
    uG=np.zeros(n);
    uT=np.zeros(n);
    
    nA=0;
    nC=0;
    nG=0;
    nT=0;

    for i in range(0,n):
        nu=seq[i];
        if nu=='A':
           uA[i]=1;
           uC[i]=0;
           uG[i]=0;
           uT[i]=0;
           nA=nA+1;
        if nu=='C':
           uC[i]=1;
           uA[i]=0;
           uG[i]=0;
           uT[i]=0;
           nC=nC+1;
        if nu=='G':
           uG[i]=1;
           uA[i]=0;
           uC[i]=0;
           uT[i]=0;
           nG=nG+1;           
        if nu=='T':
           uT[i]=1;
           uA[i]=0;
           uC[i]=0;
           uG[i]=0;
           nT=nT+1;
        
    UA=np.fft.fft(uA);
    UC=np.fft.fft(uC);
    UG=np.fft.fft(uG);
    UT=np.fft.fft(uT);
    
   
    PSA=np.abs(UA[1:])**2;
    PSC=np.abs(UC[1:])**2;
    PSG=np.abs(UG[1:])**2;
    PST=np.abs(UT[1:])**2;
    
    MA=np.zeros((1,3));
    MC=np.zeros((1,3));
    MG=np.zeros((1,3));
    MT=np.zeros((1,3));
    
    for j in range(0,3):
         MA[0,j]=nA*(n-nA)*np.sum(np.power(PSA,j+1))/(np.power(nA,j+1)*np.power(n-nA,j+1));
         MC[0,j]=nC*(n-nC)*np.sum(np.power(PSC,j+1))/(np.power(nC,j+1)*np.power(n-nC,j+1));
         MG[0,j]=nG*(n-nG)*np.sum(np.power(PSG,j+1))/(np.power(nG,j+1)*np.power(n-nG,j+1));
         MT[0,j]=nT*(n-nT)*np.sum(np.power(PST,j+1))/(np.power(nT,j+1)*np.power(n-nT,j+1));
    

    v=[MA[0,0],MC[0,0],MG[0,0],MT[0,0],MA[0,1],MC[0,1],MG[0,1],MT[0,1],MA[0,2],MC[0,2],MG[0,2],MT[0,2]]
    return v;
if __name__=='__main__':
    print(GetMomentVectorPS('ATGGAGACAGCATCACTGATAACTATACTGCTAGCAGTAACAGTAAGCAATGCAGATAAGATCTGCAT'))

