<commit_before>'''

Implementation of Lenth's Method. This returns a t-like statistic allowng the sensitivity
to be examined.

'''

import numpy as np
import matplotlib.pyplot as p
from pandas import Panel, read_hdf, DataFrame
import statsmodels.formula.api as sm

model_matrix = np.array([[1,-1,1,1,1],
                         [1,1,1,1,1],
                         [1,-1,-1,1,1],
                         [1,1,-1,1,1],
                         [1,-1,1,-1,1],
                         [1,1,1,-1,1],
                         [1,-1,-1,-1,1],
                         [1,1,-1,-1,1],
                         [1,-1,1,1,-1],
                         [1,1,1,1,-1],
                         [1,-1,-1,1,-1],
                         [1,1,-1,1,-1],
                         [1,-1,1,-1,-1],
                         [1,1,1,-1,-1],
                         [1,-1,-1,-1,-1],
                         [1,1,-1,-1,-1]
                        ])

model_matrix = DataFrame(model_matrix,columns=["Constant","Mach", "B-field", "Driving", "Temperature"])

class LenthsMethod(object):
    """This implements Lenth's Method for use in the sensitivity analysis"""
    def __init__(self, datafile, save_name=None, columns=None, rows=None, \
                 statistics=["Wavelet", "MVC", "PSpec", "Genus", "VCS", "VCA", "Tsallis", "Skewness", "Kurtosis"]):
        super(LenthsMethod, self).__init__()

        assert isinstance(datafile, str)
        self.datafile = lambda x: read_hdf(datafile,x)

        self.columns = columns
        self.rows = rows
        self.save_name = save_name
        self.statistics = statistics

        self.paramnames = None

        self.mean_respvecs = []
        self.lnvar_respvecs = []

        self.mean_fitparam = []
        self.lnvar_fitparam = []

        self.mean_signifig = []
        self.lnvar_signifig = []

    def make_response_vectors(self):

        for stat in self.statistics:
            data = self.datafile(stat)
            data = data[self.columns,self.rows]

            ## Transform to orthogonal parameterization
            data = 2*data - 1

            # Two different ways to define the responses
            # Should yield a 16 element vector for each comparison statistic
            self.mean_respvecs.append(data.mean(axis=0))
            self.lnvar_respvecs.append(np.log(data.std(axis=0)**2.))

        return self

    def fit_model(self, model=None, verbose= False):

        if model is None:
            model = "Mach*B-field*Driving*Temperature"

        for i,stat, mean, lnvar in enumerate(zip(self.statistics, self.mean_respvecs, self.lnvar_respvecs)):
            if i==0:
                model_matrix.append({"mean_resp":mean, "lnvar_resp":lnvar})
            else:
                model_matrix["mean_resp"] = mean
                model_matrix["lnvar_resp"] = lnvar

            mean_model = sm.OLS("".join(["mean_resp~",model]),data = model_matrix)
            lnvar_model = sm.OLS("".join(["lnvar_resp~",model]),data = model_matrix)

            mean_results = mean_model.fit()
            lnvar_results = lnvar_model.fit()

            self.mean_fitparam.append(mean_results.params)
            self.lnvar_fitparam.append(lnvar_results.params)

            self.paramnames = mean_model.params

            if verbose:
                print "Fits for %s" (stat)
                print mean_results.summary()
                print lnvar_results.summary()
        return self

    def lenth(self, IER=2.16):
        '''
         IER: Wu Hamada critical value 2.16 for 15 compared effect estimate with alpha 0.05
         '''

        for mean, lnvar in zip(self.mean_fitparam, self.lnvar_fitparam):
            s0_mean = np.median(np.abs(mean))
            s0_lnvar = np.median(np.abs(lnvar))

            pse_mean = 1.5*np.median(np.abs(mean[np.where(mean>2.5*s0_mean)]))
            pse_lnvar = 1.5*np.median(np.abs(lnvar[np.where(lnvar>2.5*s0_lnvar)]))

            self.mean_signifig.append([mean/pse_mean,pse_mean*IER])
            self.lnvar_signifig.append([lnvar/pse_lnvar, pse_lnvar*IER])

        return self

    def make_plots(self):
        import matplotlib.pyplot as p

        index = np.arange(len(self.statistics))
        bar_width = 0.35

        for i,stat in enumerate(self.statistics):
            p.subplot(2,1,1)
            p.bar(index, self.mean_signifig[i], color="b")
            p.xlabel("Effect Estimates")
            p.ylabel(stat)
            p.xticks(index+bar_width, self.paramnames, rotation="vertical", verticalalignment="bottom")
            p.tight_layout()
            p.subplot(2,1,2)
            p.bar(index, self.lnvar_signifig[i], color="b")
            p.xlabel("Effect Estimates")<commit_msg>Minor fixes to sensitivity analysis<commit_after>'''

Implementation of Lenth's Method. This returns a t-like statistic allowng the sensitivity
to be examined.

'''

import numpy as np
import matplotlib.pyplot as p
from pandas import Panel, read_hdf, DataFrame, concat, Series
import statsmodels.formula.api as sm


class LenthsMethod(object):
    """This implements Lenth's Method for use in the sensitivity analysis"""
    def __init__(self, datafile, save_name=None, columns=None, rows=None, \
                 statistics=["Wavelet", "MVC", "PSpec", "Genus", "VCS", "VCA", "Tsallis", "Skewness", "Kurtosis"]):
        super(LenthsMethod, self).__init__()

        assert isinstance(datafile, str)
        self.datafile = lambda x: read_hdf(datafile,x)

        self.columns = columns
        self.rows = rows
        self.save_name = save_name
        self.statistics = statistics

        self.paramnames = None

        self.mean_respvecs = []
        self.lnvar_respvecs = []

        self.mean_fitparam = []
        self.lnvar_fitparam = []

        self.mean_signifig = []
        self.lnvar_signifig = []

        self.model_matrix = np.array([[1,-1,1,1,1],
                         [1,1,1,1,1]])#,
                         # [1,-1,-1,1,1],
                        #  [1,1,-1,1,1],
                        #  [1,-1,1,-1,1],
                        #  [1,1,1,-1,1],
                        #  [1,-1,-1,-1,1],
                        #  [1,1,-1,-1,1],
                        #  [1,-1,1,1,-1],
                        #  [1,1,1,1,-1],
                        #  [1,-1,-1,1,-1],
                        #  [1,1,-1,1,-1],
                        #  [1,-1,1,-1,-1],
                        #  [1,1,1,-1,-1],
                        #  [1,-1,-1,-1,-1],
                        #  [1,1,-1,-1,-1]
                        # ])

        self.model_matrix = DataFrame(self.model_matrix,columns=["Constant","Mach", "B_field", "Driving", "Temperature"])

    def make_response_vectors(self):

        for i, stat in enumerate(self.statistics):
            data = self.datafile(stat)

            if i==0:
                self.model_matrix.index = data.index # Set the indexes to be the same

            if self.columns is not None:
                data = data[self.columns,:]
            if self.rows is not None:
                data = data[:,self.rows]

            ## Transform to orthogonal parameterization
            data = 2*data - 1

            # Two different ways to define the responses
            # Should yield a 16 element vector for each comparison statistic
            self.mean_respvecs.append(data.mean(axis=1))
            self.lnvar_respvecs.append(np.log(data.std(axis=1)**2.))

        return self

    def fit_model(self, model=None, verbose= False):

        if model is None:
            model = "Mach*B_field*Driving*Temperature"

        for i,(stat, mean, lnvar) in enumerate(zip(self.statistics, self.mean_respvecs, self.lnvar_respvecs)):
            # if i==0:
            #     # model_matrix.append({"mean_resp":mean, "lnvar_resp":lnvar})
            #     self.model_matrix = concat(self.model_matrix[:],{"mean_resp":mean, "lnvar_resp":lnvar})
            # else:
            self.model_matrix["mean_resp"] = Series(mean, index=self.model_matrix.index)
            self.model_matrix["lnvar_resp"] = Series(lnvar, index=self.model_matrix.index)

            mean_model = sm.ols("".join(["mean_resp~",model]), data=self.model_matrix)
            lnvar_model = sm.ols("".join(["lnvar_resp~",model]), data=self.model_matrix)

            mean_results = mean_model.fit()
            lnvar_results = lnvar_model.fit()

            self.mean_fitparam.append(mean_results.params[1:])
            self.lnvar_fitparam.append(lnvar_results.params[1:])

            if i==0:
                self.paramnames = mean_model.exog_names[1:] # Set the names of the coefficients

            if verbose:
                print "Fits for "+ stat
                print mean_results.summary()
                print lnvar_results.summary()
        return self

    def lenth(self, IER=2.16):
        '''
         IER: Wu Hamada critical value 2.16 for 15 compared effect estimate with alpha 0.05
         '''

        for mean, lnvar in zip(self.mean_fitparam, self.lnvar_fitparam):
            s0_mean = np.median(np.abs(mean))
            s0_lnvar = np.median(np.abs(lnvar))

            pse_mean = 1.5*np.median(np.abs(mean)[np.where(np.abs(mean)<=2.5*s0_mean)[0]])
            pse_lnvar = 1.5*np.median(np.abs(lnvar)[np.where(lnvar<=2.5*s0_lnvar)[0]])

            self.mean_signifig.append([mean/pse_mean,pse_mean*IER])
            self.lnvar_signifig.append([lnvar/pse_lnvar, pse_lnvar*IER])

        return self

    def make_plots(self):
        import matplotlib.pyplot as p

        index = np.arange(len(self.paramnames))
        bar_width = 0.35

        for i,stat in enumerate(self.statistics):
            p.subplot(2,1,1)
            p.bar(index, self.mean_signifig[i][0], color="b")
            p.plot([0.,index.max()+2*bar_width],[self.mean_signifig[i][1],self.mean_signifig[i][1]], "b-",
                   [0.,index.max()+2*bar_width],[-self.mean_signifig[i][1],-self.mean_signifig[i][1]], "b-")
            p.xlabel("Effect Estimates")
            p.ylabel(stat)
            p.xticks(index+bar_width, self.paramnames, rotation="vertical", verticalalignment="bottom")
            p.tight_layout()
            p.subplot(2,1,2)
            p.bar(index, self.lnvar_signifig[i][0], color="b")
            p.plot([0.,index.max()+2*bar_width],[self.lnvar_signifig[i][1],self.lnvar_signifig[i][1]], "b-",
                   [0.,index.max()+2*bar_width],[-self.lnvar_signifig[i][1],-self.lnvar_signifig[i][1]], "b-")
            p.xlabel("Effect Estimates")
{'fraction_non_alphanumeric': 0.14038938728765032, 'fraction_numerical': 0.02529108608513075, 'mean_word_length': 1.9165043139437796, 'pattern_counts': {'":': 6, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1802745', 'n_tokens_mistral': 3546, 'n_tokens_neox': 3350, 'n_words': 587}