"""
Function of this code: the total procedures of feature extraction from the numerical values.
Contact for help: gch_liu@163.com 
"""

from dna2numVOSS import dna2numVOSS
from dna2numInteger import dna2numInteger
from dna2numSGI import dna2numSGI
from dna2numReal import dna2numReal
from dna2numComplex import dna2numComplex
from dna2numKTPC import dna2numKTPC
from dna2numKBPC1 import dna2numKBPC1
from dna2numKQC1 import dna2numKQC1
from dna2numEIIP import dna2numEIIP
from dna2numPNAN import dna2numPNAN
from dna2numAN import dna2numAN
from dna2numMM import dna2numMM
from dna2numPN import dna2numPN
from dna2numSP import dna2numSP
from dna2numDNAwalk import dna2numDNAwalk
from dna2numHcurve import dna2numHcurve
from dna2numCGRcurve import dna2numCGRcurve
from dna2numDVcurve import dna2numDVcurve
from dna2numZcurve import dna2numZcurve
class trans:
    def __init__(self,num=[],mean=[],std=0,kurtosis=[],skewness=[],Hurst=[],fftpowerspectra=[],meanpowersp=[]):
        self.num=num;
        self.mean=mean;
        self.std=std;
        self.kurtosis=kurtosis;
        self.skewness=skewness;
        self.Hurst=[];
        self.fftpowerspectra=fftpowerspectra;
        self.meanpowersp=meanpowersp;
def feature_extract(seq_char):
    import numpy as np
    # fixed mapping techniques
    s1_VOSS=np.array(dna2numVOSS(seq_char)).T;  #VOSS
    s2_integer=np.array(dna2numInteger(seq_char)).T; #Integer Number
    s3_SGI=np.array(dna2numSGI(seq_char)).T;  #Single Galois Indicator
    s4_Real=np.array(dna2numReal(seq_char)).T;  # Real number
    s5_Complex=np.array(dna2numComplex(seq_char)).T; #Complex nmuber
    s6_KTPC=np.array(dna2numKTPC(seq_char)).T; # K-Twin-Pair Code
    s7_KBPC1=np.array(dna2numKBPC1(seq_char)).T; #K-Bipolar-Pair Code 1
    s8_KQC1=np.array(dna2numKQC1(seq_char)).T; #K-Quaternary Code 1
    
    # Physic-chemical property
    s9_EIIP=np.array(dna2numEIIP(seq_char)).T;  # the Electron-Ion Interaction Potential (EIIP)
    s10_PNAN=np.array(dna2numPNAN(seq_char)).T;  # Paired Nucleotide Atomic Number
    s11_AN=np.array(dna2numAN(seq_char)).T;  # Atomic Number
    s12_MM=np.array(dna2numMM(seq_char)).T; # Molecular Mass
    s13_PN=np.array(dna2numPN(seq_char)).T;  # Paired Numeric
    s14_SP=np.array(dna2numSP(seq_char)).T; #  Structure profile  (Florquin K, Saeys Y, Degroeve S, Rouze P, Peer YV (2005) Large-scale structural analysis of the core promoter in mammalian and plant genomes. Nucleic Acids Res 13:4255-4264.)
    s15_DNAwalk=np.array(dna2numDNAwalk(seq_char)).T; # DNA walk
    
    # DNA-graph based long range correlation methods
    s16_Hcurve=np.array(dna2numHcurve(seq_char)).T;  # H-curve  https://books.google.com/books?id=np.array(Atou18eg7E0C&pg=np.array(PA67&lpg=np.array(PA67&dq=np.array(DNA+H+curve&source=np.array(bl&ots=np.array(djyE6b0QWo&sig=np.array(g40WuHz8v3U5E-UmEtCYD6waXQc&hl=np.array(en&sa=np.array(X&ei=np.array(pgU7VfLtIIT7sASJzoCgDw&ved=np.array(0CFoQ6AEwDg#v=np.array(onepage&q=np.array(DNA#20H#20curve&f=np.array(false
    s17_CGRcurve=np.array((dna2numCGRcurve(seq_char))); # CGR curve
    s18_DVcurve=np.array((dna2numDVcurve(seq_char))).T; # Dual-Vector(DV) curve
    s19_Zcurve=np.array(dna2numZcurve(seq_char)).T; # Z curve
    # note  DVcurve can be looked as one dimension
    
    # Genomic signal processing
    trans_seq=[]
    trans_seq.append(trans(num=s1_VOSS));  
    trans_seq.append(trans(num=s2_integer));
    trans_seq.append(trans(num=s3_SGI));
    trans_seq.append(trans(num=s4_Real));
    trans_seq.append(trans(num=s5_Complex));
    trans_seq.append(trans(num=s6_KTPC));
    trans_seq.append(trans(num=s7_KBPC1));
    trans_seq.append(trans(num=s8_KQC1));
    trans_seq.append(trans(num=s9_EIIP));
    trans_seq.append(trans(num=s10_PNAN));
    trans_seq.append(trans(num=s11_AN));
    trans_seq.append(trans(num=s12_MM));
    trans_seq.append(trans(num=s13_PN));
    trans_seq.append(trans(num=s14_SP));
    trans_seq.append(trans(num=s15_DNAwalk));
    # Basic features of graph-based methods
    trans_seq.append(trans(num=s16_Hcurve));
    trans_seq.append(trans(num=s17_CGRcurve));
    trans_seq.append(trans(num=s18_DVcurve));
    trans_seq.append(trans(num=s19_Zcurve));
    
    real_trans=list(range(0,4)); # real transformation
    real_trans.extend([x for x in range(8,14)])
    real_trans.extend([x for x in range(16,19)])
    complex_trans=list(range(4,8)); # complex transformation
    complex_trans.extend([14,15])
    r1=len(real_trans);
    # c1=len(complex_trans);
    
    from GetMomentVectorPS import GetMomentVectorPS
    VOSS_moments=GetMomentVectorPS(seq_char);
    import numpy as np
    from scipy import stats
    real_feature=[]
    for j in range(0,r1):
        i=real_trans[j];
        a=trans_seq[i].mean=np.mean(trans_seq[i].num,0);
        b=trans_seq[i].std=np.std(trans_seq[i].num,0,ddof=1);
        c=trans_seq[i].kurtosis=stats.kurtosis(trans_seq[i].num,0,fisher=False);
        d=trans_seq[i].skewness=stats.skew(trans_seq[i].num,0);
        if i==17:
            a=np.array([a[1]])
            b=np.array([b[1]])
            c=np.array([c[1]])
            d=np.array([d[1]])
        if type(a)!=np.ndarray:
            a=np.array([a])
            b=np.array([b])
            c=np.array([c])
            d=np.array([d])
#        print(i,j,b)
        real_feature+=a.tolist()
        real_feature+=b.tolist()
        real_feature+=c.tolist()
        real_feature+=d.tolist()
    
           
    Hurst_exp=[]
    import genhurst
    for i in range(0,19):
        kt=(trans_seq[i].num)[0].size;
        if kt>=2:
            for j in range(0,kt):
                temp=trans_seq[i].num;
             
                a=temp[:,j]
                temphurst=genhurst.genhurst(a)
                trans_seq[i].Hurst.append(temphurst);
                if i==0 or (i==17 and j!=1):
                    pass
                else:
                    Hurst_exp.append(temphurst)

        else:
            trans_seq[i].Hurst=genhurst.genhurst(trans_seq[i].num);
            Hurst_exp.append(trans_seq[i].Hurst)

        
#% based on the FFT and power spectra, then calculate their mean,
    from seq_fft import seq_fft
    complex_feature=[]
    for i in range(0,19):
        trans_seq[i].fftpowerspectra=[]
        trans_seq[i].meanpowersp=[]
        kt=(trans_seq[i].num)[0].size;
        if kt>=2:
            for j in range(0,kt):                       
                a=np.array(trans_seq[i].num)[:,j]
                b=seq_fft(a.tolist())
                trans_seq[i].fftpowerspectra.append(b);
                temp=np.mean(trans_seq[i].fftpowerspectra[:][j])
                trans_seq[i].meanpowersp.append(temp);

                if i==17 and j==0:
                    pass
                else:
                    complex_feature.append(temp)

        else:
            trans_seq[i].fftpowerspectra=seq_fft(trans_seq[i].num);  
            temp=np.mean(trans_seq[i].fftpowerspectra);
            trans_seq[i].meanpowersp=temp;
            complex_feature.append(temp)

                     
    feature2=VOSS_moments;
    feature2+=real_feature;
    feature2+=Hurst_exp;
    feature2+=complex_feature;
    return feature2
if __name__=='__main__':
    print(feature_extract('ACGTACGTACGTACGT'))
