Ước lượng phi tham số cho P(X<Y) - Luận án tiến sĩ Tôn Thất Quang Nguyên

Bài toán xuất phát từ việc quan sát Xi = Xj + εj và Yi = Yk + ηk thay vì các giá trị thực Xj và Yk. Các biến ngẫu nhiên εj và ηk đại diện cho sai số đo. Giả định các biến Xj độc lập cùng phân phối với hàm mật độ fx chưa biết. Tương tự, các biến Yk độc lập cùng phân phối với hàm mật độ fy chưa biết. Sai số đo có thể đồng phân phối hoặc dị phân phối tùy từng trường hợp. Mục tiêu là ước lượng θ = P(X < Y) từ dữ liệu nhiễu này. Độ chính xác của ước lượng phụ thuộc vào đặc tính của sai số và kích thước mẫu.

Trường hợp thứ nhất xét sai số đồng phân phối. Các εj độc lập cùng hàm mật độ fε. Các ηk độc lập cùng hàm mật độ fη. Giá của các hàm mật độ này nằm trong [-M, M] với M > 0. Trường hợp thứ hai xét sai số dị phân phối. Mỗi εj có hàm mật độ riêng fεj. Mỗi ηk có hàm mật độ riêng fηk. Giá vẫn bị chặn trong [-M, M] với M độc lập với kích thước mẫu. Cả hai trường hợp đều yêu cầu phương pháp deconvolution để loại bỏ ảnh hưởng của sai số.

Xác suất P(X < Y) đo lường độ tin cậy hệ thống trong mô hình stress-strength. X đại diện cho độ bền của hệ thống. Y đại diện cho ứng suất tác động lên hệ thống. Hệ thống hoạt động tốt khi X > Y. Ước lượng chính xác P(X < Y) giúp đánh giá độ tin cậy. Trong thực tế, các phép đo luôn có sai số. Phương pháp phi tham số không yêu cầu giả định phân phối cứng nhắc. Điều này tăng tính linh hoạt và khả năng ứng dụng rộng rãi.

Ước lượng hàm mật độ kernel có dạng tổng các hàm nhân. Mỗi quan sát đóng góp một hàm nhân được làm trơn. Tham số bandwidth h điều chỉnh độ rộng của hàm nhân. Bandwidth nhỏ cho ước lượng chi tiết nhưng nhiễu. Bandwidth lớn cho ước lượng trơn nhưng mất thông tin. Hàm nhân thường dùng gồm Gaussian, Epanechnikov, và uniform. Tính chất của hàm nhân ảnh hưởng đến hiệu quả ước lượng. Phương pháp này không yêu cầu giả định tham số về phân phối.

Biến đổi Fourier chuyển hàm sang miền tần số. Tích chập trong miền thời gian trở thành phép nhân trong miền tần số. Nếu Z = X + ε thì φZ(t) = φX(t)φε(t) với φ là hàm đặc trưng. Ước lượng φX(t) = φZ(t)/φε(t) khi biết φε(t). Biến đổi Fourier ngược cho ước lượng hàm mật độ fx. Quá trình này gọi là deconvolution. Vấn đề ill-posed xuất hiện khi φε(t) gần 0. Cần kỹ thuật regularization để ổn định ước lượng.

Hàm mật độ trơn thường có hàm đặc trưng giảm theo lũy thừa đa thức. Hàm mật độ siêu trơn có hàm đặc trưng giảm theo hàm mũ. Sai số Laplace và Gamma là ví dụ điển hình của siêu trơn. Trường hợp siêu trơn khó hơn do deconvolution kém ổn định. Tốc độ hội tụ chậm hơn với sai số siêu trơn. Cần điều chỉnh bandwidth theo logarithm của kích thước mẫu. Phương pháp minimax giúp tìm tốc độ hội tụ tối ưu. Chặn trên và chặn dưới xác định giới hạn lý thuyết.

Bước đầu ước lượng hàm mật độ fx từ dữ liệu nhiễu Xi. Sử dụng kernel deconvolution với bandwidth hn phù hợp. Tương tự ước lượng hàm mật độ fy từ dữ liệu Yi. Bandwidth có thể khác nhau cho hai ước lượng. Tính θ̂ = ∫∫ I(x<y) f̂x(x) f̂y(y) dx dy. Tích phân này được xấp xỉ bằng phương pháp số. Sử dụng lưới điểm rời rạc để tính tổng. Độ mịn của lưới ảnh hưởng đến sai số xấp xỉ.

Tính vững yêu cầu E|θ̂n - θ|² → 0 khi n, m → ∞. Điều kiện giá compact của sai số đảm bảo khả nghịch của deconvolution. Hàm đặc trưng của sai số không triệt tiêu trên R. Bandwidth phải thỏa mãn hn → 0 và nhn → ∞. Tốc độ hội tụ của bandwidth phụ thuộc vào độ trơn của fx và fy. Với hàm mật độ thuộc lớp Sobolev, bandwidth tối ưu có dạng n^(-1/(2β+1)). Tham số β đo độ trơn của hàm mật độ. Chứng minh tính vững sử dụng bất đẳng thức Minkowski và Cauchy-Schwarz.

Ước lượng Wilcoxon-Mann-Whitney là U-statistic cổ điển. Công thức θ̂WMW = (1/nm) Σ I(Xi < Yj). Ước lượng này vững khi không có sai số đo. Với dữ liệu nhiễu, θ̂WMW ước lượng P(X+ε < Y+η). Giá trị này khác với P(X < Y) khi sai số không đối xứng. Sai lệch phụ thuộc vào phân phối của ε và η. Phương pháp deconvolution khắc phục vấn đề này. Ước lượng mới có bias nhỏ hơn nhưng variance lớn hơn. Trade-off này được kiểm soát qua lựa chọn bandwidth.

Rủi ro tối đa là sup E|θ̂n - θ|² trên lớp hàm cho trước. Phân tích sai số thành bias và variance. Bias xuất phát từ việc làm trơn qua bandwidth. Variance xuất phát từ biến động mẫu ngẫu nhiên. Với lớp Sobolev β-trơn, bias = O(h^β). Variance = O((nh)^(-1)) cho kernel bậc hai. Cân bằng bias-variance cho bandwidth tối ưu h* ~ n^(-1/(2β+1)). Thay vào được rủi ro tối đa O(n^(-2β/(2β+1))). Chặn này đúng cho cả hai ước lượng hàm mật độ.

Chặn dưới minimax chứng minh không có ước lượng nào tốt hơn. Kỹ thuật Assouad lemma xây dựng họ phân phối khó phân biệt. Chọn hai hàm mật độ fx0 và fx1 gần nhau. Khoảng cách Kullback-Leibler giữa chúng nhỏ. Nhưng θ0 = P(X0 < Y) và θ1 = P(X1 < Y) khác nhau đáng kể. Mọi ước lượng phải nhầm lẫn với xác suất không đổi. Điều này dẫn đến rủi ro tối thiểu không thể tránh. Chặn dưới phù hợp với chặn trên đã thiết lập. Kết luận tốc độ hội tụ là tối ưu.

Sai số Laplace có hàm mật độ fε(x) = (1/2)e^(-|x|). Hàm đặc trưng φε(t) = 1/(1+t²) giảm chậm. Deconvolution yêu cầu chia cho φε(t) rất nhỏ ở tần số cao. Điều này khuếch đại nhiễu nghiêm trọng. Bandwidth phải chọn hn ~ (log n)^(-1) thay vì n^(-1/(2β+1)). Tốc độ hội tụ trở thành (log n)^(-β) chậm hơn nhiều. Sai số Gamma có tính chất tương tự. Kết quả cho thấy sai số siêu trơn khó xử lý hơn đáng kể. Tuy nhiên ước lượng vẫn hội tụ về giá trị thực.

Bước 1: Từ {X1,...,Xn} rút mẫu có hoàn lại {X1,...,Xn}. Bước 2: Từ {Y1,...,Ym} rút mẫu có hoàn lại {Y1,...,Ym}. Bước 3: Tính ước lượng θ̂b từ mẫu bootstrap. Sử dụng cùng bandwidth và phương pháp deconvolution. Bước 4: Lặp lại B lần để có {θ̂1,...,θ̂B}. Thông thường B = 500 hoặc 1000 là đủ. Bước 5: Tính phương sai mẫu của các θ̂b. Bước 6: Tính quantile 2.5% và 97.5% cho khoảng tin cậy 95%.

Bootstrap nhất quán nếu phân phối bootstrap hội tụ đến phân phối giới hạn. Với dữ liệu không nhiễu, điều kiện thông thường đủ. Với dữ liệu nhiễu, cần điều kiện bổ sung về sai số. Sai số phải có moment bậc cao hữu hạn. Bandwidth bootstrap phải giống bandwidth gốc. Tốc độ hội tụ của bootstrap có thể chậm hơn ước lượng gốc. Trong một số trường hợp, cần m-out-of-n bootstrap. Phương pháp này lấy mẫu kích thước nhỏ hơn n. Điều này cải thiện tính chất của bootstrap với dữ liệu nhiễu.

Kiểm định H0: θ = θ0 với mức ý nghĩa α. Tính khoảng tin cậy (1-α) cho θ. Bác bỏ H0 nếu θ0 nằm ngoài khoảng tin cậy. Bootstrap cung cấp p-value xấp xỉ. P-value = tỷ lệ |θ̂*b - θ̂n| > |θ0 - θ̂n|. So sánh hai hệ thống qua H0: θ = 0.5. θ > 0.5 nghĩa hệ thống X tin cậy hơn Y. θ < 0.5 nghĩa hệ thống Y tin cậy hơn X. Phương pháp bootstrap linh hoạt cho nhiều dạng kiểm định.

Chọn fx là phân phối chuẩn N(μx, σx²). Chọn fy là phân phối chuẩn N(μy, σy²). Điều chỉnh μx, μy để thay đổi θ thực. Sai số ε ~ N(0, σε²) cho trường hợp trơn thường. Sai số η ~ Laplace(0, b) cho trường hợp siêu trơn. Kích thước mẫu n = m thay đổi từ 100 đến 1000. Số lần lặp Monte Carlo là 1000. Tính bias = trung bình(θ̂ - θ). Tính MSE = trung bình(θ̂ - θ)². So sánh với ước lượng Wilcoxon-Mann-Whitney.

Ước lượng deconvolution có bias nhỏ hơn đáng kể so với WMW. Với sai số không đối xứng, WMW sai lệch lên đến 0.1. Phương pháp đề xuất giữ bias dưới 0.01. MSE giảm theo tốc độ lý thuyết khi n tăng. Với sai số trơn thường, MSE ~ n^(-0.8) phù hợp lý thuyết. Với sai số siêu trơn, MSE giảm chậm hơn. Lựa chọn bandwidth ảnh hưởng lớn đến kết quả. Cross-validation cho bandwidth gần tối ưu. Độ phủ của khoảng tin cậy bootstrap đạt 94-96%. Thời gian tính toán chấp nhận được với n < 500.

Xét hệ thống điện tử với tuổi thọ X. Ứng suất hoạt động trong thời gian Y. Đo lường tuổi thọ có sai số do điều kiện thử nghiệm. Đo lường ứng suất có sai số do cảm biến. Dữ liệu gồm 200 quan sát tuổi thọ và 150 quan sát ứng suất. Ước lượng WMW cho θ̂ = 0.73. Phương pháp deconvolution cho θ̂ = 0.68. Khoảng tin cậy 95% là [0.62, 0.74]. Kết luận hệ thống có độ tin cậy hệ thống khoảng 68%. Điều chỉnh cho sai số thay đổi kết luận đáng kể. Quyết định thiết kế dựa trên ước lượng chính xác hơn.