Değişkenlik hipotezi#

Değişkenlik hipotezi, erkeklerin fiziksel ve psikolojik özelliklerde genel olarak kadınlardan daha değişken olduğu kuramıdır. Hafif söylemek gerekirse tartışmalı bir iddiadır.

Konuyu hızlıca incelemek ve bootstrap yeniden örneklemesini göstermek için erkeklerin boylarının kadınlardan daha değişken olup olmadığına bakalım. 400.000’den fazla katılımcının kendi bildirdiği boy ve ağırlıkları içeren 2022 BRFSS verisini kullanacağız.

try:
    import empiricaldist
except ImportError:
    !pip install empiricaldist
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

np.random.seed(42)
plt.rcParams["figure.dpi"] = 75
plt.rcParams["figure.figsize"] = [6, 3.5]

Veriler#

Aşağıdaki işlev, 1 Ağustos 2024’te BRFSS sitesinden indirilmiş SAS dışa aktarma biçimindeki veriyi okur. Kullanacağımız sütunları seçip çok daha küçük alt kümeyi HDF dosyasına yazar.

def write_hdf():
    """SAS dışa aktarma dosyasını okuyup bir HDF dosyası yazar."""
    brfss = pd.read_sas("LLCP2022.XPT.gz")
    columns = ["_SEX", "HTM4", "WTKG3", "_LLCPWT"]
    brfss[columns].to_hdf("LLCP2022.hdf5", key="brfss", complevel=6)

Tüm veri kümesi yerine yalnızca alt kümeyi indirebiliriz.

from os.path import basename, exists


def download(url):
    filename = basename(url)
    if not exists(filename):
        from urllib.request import urlretrieve

        local, _ = urlretrieve(url, filename)
        print("İndirildi: " + local)


download("https://github.com/TerekliTahaBerk/thinkstatstr/raw/v3/data/LLCP2022.hdf5")

İlk birkaç satır şöyledir.

brfss = pd.read_hdf("LLCP2022.hdf5", key="brfss")
brfss.describe()
_SEX HTM4 WTKG3 _LLCPWT
count 445132.000000 416480.000000 403054.000000 445132.000000
mean 1.529942 170.269057 8307.447039 594.856344
std 0.499103 10.717750 2144.817270 1134.837415
min 1.000000 91.000000 2268.000000 0.020464
25% 1.000000 163.000000 6804.000000 115.885991
50% 2.000000 170.000000 8074.000000 274.632388
75% 2.000000 178.000000 9525.000000 627.913694
max 2.000000 241.000000 29257.000000 54390.520926

Veri kümesi yaklaşık 200.000 erkek ve 220.000 kadının bildirdiği boyları içerir. Kaydedilen en kısa boylara sahip erkek ve kadın sayıları şöyledir.

xtab = pd.crosstab(brfss["HTM4"], brfss["_SEX"])
xtab.columns = ["Erkek", "Kadın"]
xtab.head()
Erkek Kadın
HTM4
91.0 7 17
92.0 0 1
95.0 1 0
97.0 1 3
99.0 1 0

En uzun kaydedilen boylar da şöyledir.

xtab.tail()
Erkek Kadın
HTM4
226.0 9 2
229.0 4 1
234.0 4 0
236.0 1 0
241.0 4 1

Kod kitabına göre 91 cm altındaki boylar 91, 241 cm üzerindekiler 241 cm yapılır. Yine de kalan bazı değerler muhtemelen hatalıdır. Dünyanın yaşayan en uzun kadını 234 cm olduğundan veri kümesindeki 241 cm büyük olasılıkla yanlış bildirilmiş ya da kaydedilmiştir.

Bu DataFrame’den bootstrap örneklemi çekmek için aşağıdaki işlevi kullanacağız.

def resample(df):
    """Bir bootstrap örneklemi çeker.

    df: DataFrame

    döndürür: DataFrame
    """
    n = len(df)
    return df.sample(n, replace=True, weights="_LLCPWT")

Sınamada kullanacağımız tek bir örneklem şöyledir.

sample = resample(brfss)

Ortalamalar arasındaki fark#

Erkeklerin daha değişken olup olmadığını sınamadan önce ortalama olarak daha uzun olduklarını doğrulayalım. Ortalama farkının örnekleme dağılımını bootstrap ile tahmin edeceğiz. Aşağıdaki yardımcı işlev bir dizinin ilk iki öğesi arasındaki farkı hesaplar.

def diff(seq):
    """Bir dizinin ilk iki öğesi arasındaki fark.

    seq: dizi

    döndürür: float
    """
    return np.diff(seq)[0]

Bunu bir örneklemdeki ortalama farkını, yani “test istatistiğini” hesaplamak için kullanacağız.

def diff_means(sample):
    """Ortalama boy farkı (erkek eksi kadın).

    sample: _SEX ve HTM4 sütunlarını içeren DataFrame

    döndürür: float
    """
    grouped_heights = sample.groupby("_SEX")["HTM4"]
    return -diff(grouped_heights.mean())

Fark, pozitif değerin erkeklerin daha uzun olduğunu göstereceği biçimde tanımlıdır. Bootstrap örneklemiyle örnek şöyledir.

diff_means(sample)
np.float64(14.455796580284954)

Erkekler ortalama yaklaşık 14 cm daha uzundur.

Aşağıdaki işlev bir DataFrame ve test istatistiği hesaplayan işlev alır. Çok kez yeniden örnekler, her örneklemin istatistiğini hesaplar ve liste döndürür.

def sampling_dist(df, test_stat, iters=201):
    """Bootstrap örneklemleri üretip test istatistiğini hesaplar.

    df: DataFrame
    test_stat: DataFrame alan işlev
    iters: örneklem sayısı

    döndürür: float listesi
    """
    return [test_stat(resample(df)) for i in range(iters)]

Şöyle çağırırız. Sonuç ortalama farklarının örnekleme dağılımından bir örneklemdir.

t1 = sampling_dist(brfss, diff_means)
np.mean(t1)
np.float64(14.433755828000592)

Aşağıdaki işlev KDE ile örnekleme dağılımının yoğunluğunu tahmin eder, güven aralığını hesaplar ve dağılımı, ortalamayı ve %95 GA’yı gösteren grafik üretir.

from scipy.stats import gaussian_kde
from empiricaldist import Pmf


def plot_sampling_distribution(data, format_str=".2f", **options):
    """Örnekleme dağılımını, ortalamayı ve güven aralığını çizer.

    data: örnekleme dağılımından örneklem
    format_str: değerleri biçimlendirmek için kullanılır
    options: plt.plot'a iletilen seçenekler
    """
    # ortalamayı ve güven aralığını hesapla
    m, s = np.mean(data), np.std(data)
    qs = np.linspace(m - 4 * s, m + 4 * s, 501)
    ps = gaussian_kde(data)(qs)
    pmf = Pmf(ps, qs)
    pmf.normalize()
    low, high = pmf.make_cdf().inverse([0.025, 0.975])

    # yoğunluğu çiz
    pmf.plot(**options)
    inside = pmf[(pmf.qs >= low) & (pmf.qs <= high)]
    plt.fill_between(inside.index, inside, color="gray", alpha=0.2)

    y = pmf.max()
    ax = plt.gca()

    # metinleri ve çizgileri ekle
    plt.text(m, 0.4 * y, "ortalama", fontsize=14, ha="center")
    plt.text(m, 0.05 * y, f"%95 GA", fontsize=14, ha="center")
    plt.plot([m, m], [0, y], ls=":", color="gray")
    plt.plot([low, high], [0, 0], ls="-", color="gray")

    # adjust xticks
    ticks = [low, m, high]
    labels = [format(tick, format_str) for tick in ticks]
    ax.tick_params(axis="x", labelsize=14)
    plt.xticks(ticks, labels)
    plt.yticks([])

    # remove spines
    for spine in ax.spines.values():
        spine.set_visible(False)

Ortalama farkının örnekleme dağılımı şöyledir.

plot_sampling_distribution(t1)
plt.title("Ortalama farkı (cm)", fontsize=14)
plt.tight_layout()
plt.savefig("variability1.png", dpi=300)
_images/f0c6d6f5e4f7ddedb6b2e79016f9211ca8a9e0d97063d7ac2445bb80d8f90dcc.png

Ortalama farkı yaklaşık 14,4 cm’dir. Veri toplama sürecini çok kez simüle edince fark biraz değişir, ancak sıfıra yaklaşmaz.

Aşağıdaki işlev bu dağılımdan bir değerin 0’ı aşma olasılığını, yani erkeklerin ortalama daha uzun olduğu hipotezinin p-değerini tahmin eder.

from scipy.stats import norm


def p_value(t):
    m, s = np.mean(t), np.std(t)
    if m > 0:
        return norm.cdf(0, m, s)
    else:
        return norm.sf(0, m, s)

p-değeri kayan nokta aritmetiğiyle hesaplanamayacak kadar sıfıra yakındır; gerçekten ihmal edilebilir.

p_value(t1)
np.float64(0.0)

Aynı evrenden bu büyüklükte başka bir örneklem alırsak kadınların ortalama olarak daha uzun çıkması uygulamada olanaksızdır.

Örneklemde gördüğümüz farkın şans eseri oluşmasının çok düşük olasılıklı olduğu sonucuna da makul biçimde varabiliriz. Hipotez testlerini çok biçimsel yorumlayan bazı kişiler bunu sevmez, ancak görüşümün arkasındayım.

Standart sapmalar arasındaki fark#

Şimdi asıl soruya gelelim. Erkekler daha değişkense boylarının standart sapmasının daha yüksek olmasını bekleriz. Bir örneklem alıp standart sapma farkını hesaplayan işlev yeterlidir.

def diff_stds(sample):
    """Standart sapma farkı (erkek eksi kadın).

    sample: `_SEX` ve `HTM4` sütunlarını içeren DataFrame

    döndürür: float
    """
    grouped_heights = sample.groupby("_SEX")["HTM4"]
    return -diff(grouped_heights.std())

Tek yeniden örneklemeye göre erkek boylarının standart sapması daha yüksek görünüyor.

diff_stds(sample)
np.float64(0.7755433526023632)

sampling_dist işlevini farklı test istatistiğiyle yeniden kullanabiliriz.

t2 = sampling_dist(brfss, diff_stds)
np.mean(t2)
np.float64(0.8125943377103155)

Ortalaması ve %95 GA’sıyla örnekleme dağılımı şöyledir.

plot_sampling_distribution(t2)
plt.title("Standart sapma farkı (cm)", fontsize=14)
plt.tight_layout()
plt.savefig("variability2.png", dpi=300)
_images/baf6230138f64e51b2689b1eb7f1ec00692f2af631a5e99a1ea0242efc61ebd3.png
p_value(t2)
np.float64(3.3148844524687797e-169)

p-değeri ancak hesaplanabilecek kadar büyüktür; yine de gözlenen farkın şanstan kaynaklanmasının çok düşük olasılıklı olduğu sonucuna varacak kadar küçüktür.

Ancak standart sapmaları karşılaştırmak erkeklerin daha değişken olup olmadığını değerlendirmek için anlamlı olmayabilir. Erkekler ortalama daha uzundur; mutlak olarak daha fazla değişmeleri göreli olarak da daha fazla değiştikleri anlamına gelmez. Ortalamaya göre değişim daha ilgiliyse değişim katsayısını ele alabiliriz.

Değişim katsayıları arasındaki fark#

Değişim katsayısı standart sapmanın ortalamaya oranıdır; değişkenliği göreli olarak nicelleştirir. Aşağıdaki işlev erkeklerle kadınların değişim katsayısı farkını hesaplar.

def diff_cvs(sample):
    """DK farkı (erkek eksi kadın).

    sample: `_SEX` ve `HTM4` sütunlarını içeren DataFrame

    döndürür: float
    """
    grouped_heights = sample.groupby("_SEX")["HTM4"]
    means = grouped_heights.mean()
    stds = grouped_heights.std()
    return -diff(stds / means)

Farkların örnekleme dağılımından bir örneklem şöyledir.

t3 = sampling_dist(brfss, diff_cvs)
np.mean(t3)
np.float64(0.000673958474897113)

Örnekleme dağılımı ve %95 GA şöyledir.

plot_sampling_distribution(t3, format_str="0.4f")
plt.title("Değişim katsayısı farkı", fontsize=14)
plt.tight_layout()
plt.savefig("variability3.png", dpi=300)
_images/004e68c868200087b142e9200e10acae4eaa7fef5f9ea155861b6c24cc8a8e12.png

GA sıfırı içermez; fark istatistiksel olarak anlamlıdır. p-değeri de küçüktür.

p_value(t3)
np.float64(4.91233303730651e-05)

Ancak uygulamada fark çok küçüktür. Değişim katsayıları yaklaşık 0,048’dir.

grouped_heights = sample.groupby("_SEX")["HTM4"]
means = grouped_heights.mean()
stds = grouped_heights.std()
CVs = stds / means
CVs
_SEX
1.0    0.048558
2.0    0.048107
Name: HTM4, dtype: float64

DK farkını katsayıların yüzdesi olarak yazarsak yalnızca yaklaşık %1,4’tür.

np.mean(t3) / CVs * 100
_SEX
1.0    1.387945
2.0    1.400971
Name: HTM4, dtype: float64

Bu ölçüye göre erkek boyları biraz daha değişkendir; farkın uygulamada sonucu olması beklenmez.

Ayrıca standart sapma ile değişim katsayısı aykırı değerlerden etkilenir ve veri kümesindeki uçların muhtemelen hatalı olduğunu gördük. Görünürdeki fark veri hatalarının sonucu olabilir.

Daha dirençli istatistikle yeniden deneyelim.

Dirençli DK farkı#

Değişim katsayısına alternatif olarak daha dirençli değişkenlik ölçüsü kullanabiliriz. Seçeneklerden biri medyan mutlak sapmanın (MAD) medyana oranıdır. Aşağıdaki işlev bunu hesaplar.

def robust_cv(series):
    """Medyan mutlak sapma.

    series: sayı serisi

    döndürür: float
    """
    m = series.median()
    deviations = series - m
    return deviations.abs().median() / m

Aşağıdaki işlev dirençli DK farkını hesaplar. İlk satır verilere titreşim ekler; aksi durumda her yeniden örneklemede aynı sonuçları alırdık.

def diff_robust_cv(sample):
    """Dirençli DK farkı (erkek eksi kadın).

    sample: `_SEX` ve `HTM4` sütunlarını içeren DataFrame

    döndürür: float
    """
    sample["HTM4"] += np.random.normal(0, 1, size=len(sample))
    grouped_heights = sample.groupby("_SEX")["HTM4"]
    rcvs = grouped_heights.apply(robust_cv)
    return -diff(rcvs)

Örneklemin dirençli DK farkı negatiftir; kadın boylarının daha değişken olduğunu düşündürür.

diff_robust_cv(sample)
np.float64(-0.0025971928221166715)

Örnekleme dağılımı şöyledir.

t4 = sampling_dist(brfss, diff_robust_cv)
np.mean(t4)
np.float64(-0.002502441445045549)
plot_sampling_distribution(t4, format_str="0.4f")
plt.title("Dirençli DK farkı", fontsize=14)
plt.tight_layout()
plt.savefig("variability4.png", dpi=300)
_images/7ef5c605d8f77d4759b2c2f7658b52f92e1aa84872e135a90bfb2fed72aa32c1.png

p-değeri de şöyledir.

p_value(t4)
np.float64(6.93074075784399e-122)

Dirençli DK’ye göre kadınlar daha değişkendir; fakat fark yine uygulamada önemsizdir.

Tartışma#

Bu örnek bootstrap yeniden örneklemesinin gücünü gösterir: birden çok test istatistiğini incelemek ve hipoteze en uygun olanı değerlendirmek kolaydır.

  • Boyların standart sapmasını karşılaştırırsak erkekler daha değişkendir; ama daha uzun olduklarından bu şaşırtıcı değildir. Büyük şeyler üreten süreçler çoğu zaman daha değişken şeyler de üretir.

  • Farklı büyüklükleri karşılaştırmak için değişim katsayısı daha iyi olabilir. Buna göre erkekler yine biraz daha değişkendir; ancak fark uygulamada önemsizdir.

  • Sonuç veri hatalarının ürünü olabilir. DK’nin dirençli alternatifi MAD/medyan kullanılırsa kadınlar daha değişkendir; fakat yine önemli bir fark değildir.

Sonuç olarak insan boyu dağılımı değişkenlik hipotezine güçlü destek sağlamaz.

ÇDA / medyan farkı#

Bir dirençli değişkenlik ölçüsü daha.

def diff_relative_iqr(sample):
    """ÇAA / medyan farkı (erkek eksi kadın).

    sample: `_SEX` ve `HTM4` sütunlarını içeren DataFrame

    döndürür: float
    """
    sample["HTM4"] += np.random.normal(0, 1, size=len(sample))
    grouped_heights = sample.groupby("_SEX")["HTM4"]
    medians = grouped_heights.median()
    
    # sonuçlar kullandığımız kantillere bağlıdır
    iqrs = grouped_heights.quantile(0.75) - grouped_heights.quantile(0.25)
    #iqrs = grouped_heights.quantile(0.95) - grouped_heights.quantile(0.05)
    
    return -diff(iqrs / medians)
diff_relative_iqr(sample)
np.float64(-0.0031180661772347157)
t5 = sampling_dist(brfss, diff_relative_iqr)
np.mean(t5)
np.float64(-0.004621824998497872)
plot_sampling_distribution(t5, format_str="0.4f")
plt.title("ÇAA / medyan farkı", fontsize=14)
plt.tight_layout()
plt.savefig("variability5.png", dpi=300)
_images/02b36ea3812f38f7670dc8426efc0aacbf4c0812dda9665a3b48f9cfda1e40d8.png
p_value(t5)
np.float64(6.773833573551112e-98)