報告二サンプルK-S統計

Question

問題：

ここでは、それぞれが218億個のデータポイントを含む（リストdatasetに）テキストファイルに格納された2つのデータセットをプロットします。これは、データを大き過ぎてメモリとして配列として保持することができません。私はまだそれらをヒストグラムとしてグラフにすることができますが、2 sample KS testでその差を計算する方法は不明です。これは、pltオブジェクト内の各ヒストグラムにアクセスする方法がわからないからです。

例：ここでは

はダミーデータを生成するためのいくつかのコードです：

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') 

質問：

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可能にした）を生成

ほとんどのexamples for KS-statisticsは、生データ/観測/ポイント/などの2つの配列を使用していますが、このアプローチを使用するのに十分なメモリがありません。上記の例につき、どのように私は2つのディストリビューション間のKS統計量を計算するために'gsl_test_1.txt'と'gsl_test_2.txt'から（これらの事前に計算ビンにアクセスすることができます

ボーナスカルマ：？！ 録音KS統計量とp値グラフ上

Answer 1

あなたのコードを整理しましたStringIO b/cに書き込むよりもファイルに書き込むよりも効率的です。seaborn b/c matplotlibデフォルトは不公平ですyoure bins必要に応じて両方のサンプルでしきい値を同じにする必要があります私はあなたが反復し、すべてのことをこのようにしなければならないと思いますそれは必要以上に時間がかかるでしょう。あなたは実際にそのCounterものを試してみる必要があります。一度ループすれば...同じビンサイズにすることができます。とにかくそれらを一緒にビニングしているので、浮動小数点を整数に変換するだけです。 from collections import Counterであり、C = Counter()およびC[value] += 1である。あなたはlist(C.keys())からビンを作ることができる最後にdictを持っています。あなたのデータは非常にぎこちないので、これは良いでしょう。また、pandasb/cの代わりにchunksizeをnumpyとする方法があるかどうかを確認する必要があります。numpyはインデックス作成の方が早いでしょう。 を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]))