2016-05-16 2 views
1

Проблемы:Сообщите двухвыборочную K-S статистика из двух предварительно вычисленных гистограмм

Здесь я сюжет 2 наборов данных, хранящихся в текстовых файлах (в списке dataset) каждый из которых содержит 21,8 млрд точек данных. Это делает данные слишком большими для хранения в памяти как массива. Я все еще могу нарисовать их как гистограммы, но я не уверен, как рассчитать их разницу через 2 sample KS test. Это связано с тем, что я не могу понять, как получить доступ к каждой гистограмме в объекте plt.

Пример:

Вот код, чтобы генерировать фиктивные данные:

mu = [100, 120] 
sigma = 30 
dataset = ['gsl_test_1.txt', 'gsl_test_2.txt'] 
for idx, file in enumerate(dataset): 
    dist = np.random.normal(mu[idx], sigma, 10000) 
    with open(file, 'w') as g: 
     for s in dist: 
      g.write('{}\t{}\t{}\n'.format('stuff', 'stuff', str(s))) 

Это генерирует два моих гистограмм (стало возможным here):

chunksize = 1000 
dataset = ['gsl_test_1.txt', 'gsl_test_2.txt'] 
for fh in dataset: 
    # find the min, max, line qty, for bins 
    low = np.inf 
    high = -np.inf 

    loop = 0 
    for chunk in pd.read_table(fh, header=None, chunksize=chunksize, delimiter='\t'): 
     low = np.minimum(chunk.iloc[:, 2].min(), low) 
     high = np.maximum(chunk.iloc[:, 2].max(), high) 
     loop += 1 
    lines = loop*chunksize 

    nbins = math.ceil(math.sqrt(lines)) 

    bin_edges = np.linspace(low, high, nbins + 1) 
    total = np.zeros(nbins, np.int64) # np.ndarray filled with np.uint32 zeros, CHANGED TO int64 

    for chunk in pd.read_table(fh, header=None, chunksize=chunksize, delimiter='\t'): 

     # compute bin counts over the 3rd column 
     subtotal, e = np.histogram(chunk.iloc[:, 2], bins=bin_edges) # np.ndarray filled with np.int64 

     # accumulate bin counts over chunks 
     total += subtotal 


    plt.hist(bin_edges[:-1], bins=bin_edges, weights=total) 
    plt.savefig('gsl_test_hist.svg') 

Вопрос:

Большинство examples for KS-statistics используют два массива необработанных данных/наблюдений/точек/и т. Д., Но у меня недостаточно памяти для использования этого подхода. За приведенном выше примере, как я могу получить доступ к этим предвычисленными бункеров (от 'gsl_test_1.txt' и 'gsl_test_2.txt' вычислить статистику KS между двумя распределениями

Bonus карму:?! Запись статистики KS и p-значение на графике

ответ

1

Я очистил ваш код. Запись в StringIO b/c более эффективна, чем запись в файл. Установите значение по умолчанию vibe w/seaborn b/c matplotlib defaults are haggard. Пороги youre bins должны быть одинаковыми для обоих образцов, если вы хотите тест stat для выравнивания. Я не думаю, что вы должны перебирать и делать бункеры таким образом все это будет идти дольше, чем нужно. вы действительно должны попробовать эту вещь Counter, поэтому вам нужно всего лишь прокрутить один раз ... плюс вы сможете сделать тот же размер бункера. просто конвертируйте float в ints, так как вы все равно их объединяете. from collections import Counter затем C = Counter() и C[value] += 1. у вас будет dict в конце, где вы можете сделать бункеры с list(C.keys()). это было бы хорошо, так как ваши данные настолько суровы. Кроме того, вы должны увидеть, есть ли способ сделать chunksize с numpy вместо pandas b/c numpy ПУТЬ быстрее при индексировании. попробуйте %timeit для DF.iloc[i,j] и ARRAY[i,j], и вы увидите, что я имею в виду. я сделал все, что в функции, так что это более модульным

import numpy as np 
import pandas as pd 
import math 
import matplotlib.pyplot as plt 
from io import StringIO 
from scipy.stats import ks_2samp 
import seaborn as sns; sns.set() 

%matplotlib inline 

#Added seaborn b/c it looks mo betta 

mu = [100, 120] 
sigma = 30 

def write_random(file,mu,sigma=30): 
    dist = np.random.normal(mu, sigma, 10000) 
    for i,s in enumerate(dist): 
     file.write('{}\t{}\t{}\n'.format("label_A-%d" % i, "label_B-%d" % i, str(s))) 
    return(file) 

#Writing to StringIO instead of an actual file 
gs1_test_1 = write_random(StringIO(),mu=100) 
gs1_test_2 = write_random(StringIO(),mu=120) 

chunksize = 1000 

def make_hist(fh,ax): 
    # find the min, max, line qty, for bins 
    low = np.inf 
    high = -np.inf 

    loop = 0 

    fh.seek(0) 
    for chunk in pd.read_table(fh, header=None, chunksize=chunksize, sep='\t'): 
     low = np.minimum(chunk.iloc[:, 2].min(), low) #btw, iloc is way slower than numpy array indexing 
     high = np.maximum(chunk.iloc[:, 2].max(), high) #you might wanna import and do the chunks with numpy 
     loop += 1 
    lines = loop*chunksize 

    nbins = math.ceil(math.sqrt(lines)) 

    bin_edges = np.linspace(low, high, nbins + 1) 
    total = np.zeros(nbins, np.int64) # np.ndarray filled with np.uint32 zeros, CHANGED TO int64 

    fh.seek(0) 
    for chunk in pd.read_table(fh, header=None, chunksize=chunksize, delimiter='\t'): 

     # compute bin counts over the 3rd column 
     subtotal, e = np.histogram(chunk.iloc[:, 2], bins=bin_edges) # np.ndarray filled with np.int64 

     # accumulate bin counts over chunks 
     total += subtotal 

    plt.hist(bin_edges[:-1], bins=bin_edges, weights=total,axes=ax,alpha=0.5) 

    return(ax,bin_edges,total) 

#Make the plot canvas to write on to give it to the function 
fig,ax = plt.subplots() 

test_1_data = make_hist(gs1_test_1,ax) 
test_2_data = make_hist(gs1_test_2,ax) 

#test_1_data[1] == test_2_data[1] The bins should be the same if you're going try and compare them... 
ax.set_title("ks: %f, p_in_the_v: %f" % ks_2samp(test_1_data[2], test_2_data[2])) 

enter image description here

Смежные вопросы