intTypePromotion=1
zunia.vn Tuyển sinh 2024 dành cho Gen-Z zunia.vn zunia.vn
ADSENSE

Giáo trình Thiết kế thí nghiệm: Phần 2 (Dùng cho giảng dạy cao học các ngành Chăn nuôi và Chăn nuôi – Thú y)

Chia sẻ: _ _ | Ngày: | Loại File: PDF | Số trang:83

13
lượt xem
4
download
 
  Download Vui lòng tải xuống để xem tài liệu đầy đủ

Tiếp nội dung phần 1, Giáo trình Thiết kế thí nghiệm Phần 2 do Đỗ Đức Lực và Hà Xuân Bộ biên soạn, với các nội dung chính như sau: Tương quan và hồi quy; kiểm định và so sánh tỷ lệ. Mời các bạn cùng tham khảo!

Chủ đề:
Lưu

Nội dung Text: Giáo trình Thiết kế thí nghiệm: Phần 2 (Dùng cho giảng dạy cao học các ngành Chăn nuôi và Chăn nuôi – Thú y)

  1. Chương 3. TƯƠNG QUAN VÀ HỒI QUY Đánh giá mối liên hệ của các tính trạng sản xuất trong chăn nuôi góp phần đưa ra những giải pháp phù hợp để cải thiện năng suất và nâng cao hiệu quả chăn nuôi. Nội dung của chương này đề cập đến đánh giá mối liên hệ giữa các tính trạng định lượng hoặc ảnh hưởng của một hoặc nhiều yếu tố thí nghiệm đến các chỉ tiêu nghiên cứu dưới dạng định lượng bằng phân tích hồi quy tuyến tính và phi tuyến tính. Bên cạnh đó, chương này cũng cung cấp kiến thức phương pháp bề mặt đáp ứng, đây là phương pháp nghiên cứu tối ưu hóa các yếu tố thí nghiệm trong các nghiên cứu dinh dưỡng và vi sinh vật. 3.1. TƯƠNG QUAN Tương quan được áp dụng đối với các biến định lượng và thể hiện mối liên hệ giữa 2 biến. Phân tích hồi quy được áp dụng khi cả biến phụ thuộc (số liệu thí nghiệm) có dạng biến liên tục và biến độc lập (nhân tố ảnh hưởng) cũng thuộc dạng biến liên tục. Dựa trên lý thuyết xác suất về hệ số tương quan, chúng ta có công thức sau để tính hệ số tương quan mẫu rXY giữa hai biến ngẫu nhiên X và Y: n  ( x − x )( y − y ) i i CovXY rXY = 1 = n n s X sY ( x − x ) ( y − y ) 2 2 i i 1 1 Ví dụ 3.1: Tiến hành cân khối lượng (P), đo đường kính lớn (D) và đường kính bé (d) của 22 quả trứng gà. Số liệu thu được trình bày ở bảng dưới đây. P (g) 66,80 60,10 71,20 61,60 61,20 59,00 67,90 59,00 51,50 62,60 64,20 D (mm) 58,37 54,95 60,58 56,73 57,36 53,26 57,07 58,17 52,28 55,62 56,82 d (mm) 45,12 44,35 45,56 44,34 43,57 44,86 46,27 42,82 41,91 44,95 44,79 P (g) 71,20 54,20 54,50 69,10 55,90 66,00 68,00 62,00 56,70 67,00 53,80 D (mm) 61,15 54,24 54,99 60,99 54,41 58,19 59,93 56,80 55,66 58,49 52,44 d (mm) 46,00 42,58 42,32 44,85 42,62 45,69 45,50 44,20 42,41 45,56 43,38 Cấu trúc số liệu Số liệu được nhập vào trong cửa sổ Worksheet của phần mềm Minitab; bao gồm 3 cột: 111
  2. 1) Cột Khối lượng trứng C1 (P); 2) Cột Đường kính lớn C2 (DKL); 3) Cột Đường kín bé C3 (DKN). Tính hệ số tương quan giữa Khối lượng với đường kính lớn và đường kính nhỏ được thực hiện như sau: Stat ➔ Basic Statistics ➔ Correlation... Khai báo các biến cần tính hệ số tương quan P, DKL và DKN vào ô Variables. Chọn OK để có kết quả. Kết quả out put trong cửa sổ Session của Minitab: 112
  3. Correlations: P, DKL, DKN  P DKL DKL 0.897 0.000 DKN 0.905 0.648 0.000 0.001 Cell Contents: Pearson correlation P-Value Trong đó:  Kết quả hệ số tương quan giữa 3 biến (P, DKL và DKN);  Giá trị phía trên là hệ số tương quan (r = 0,897);  Giá trị phía dưới là xác suất đối với hệ số tương quan (P < 0,000). 3.2. HỒI QUY 3.2.1. Hồi quy tuyến tính đơn và đa biến Có thể xây dựng phương trình hồi quy đơn biến y = a + bx hoặc đa biến y = a + b1x1 + b2x2 + ... +bnxn. Với ví dụ sau, ta có thể xây dựng phương trình hồi tuyến tính đơn biến quy ước tính khối lượng trứng thông qua đường kính lớn/đường kính bé hoặc đa biến thông qua đường kính lớn và đường kính bé. Ví dụ 3.1 được sử dụng lại để xây dựng phương trình hồi quy tuyến tính. Xây dựng phương trình hồi quy tuyến tính với biến phụ thuộc Y là khối lượng trứng và biến độc lập X là đường kính lớn. Mô hình phân tích: Y = a + bX +  Trong đó: + Y: Biến phụ thuộc - khối lượng trứng (g); + a: Tung độ gốc (điểm chắn); + b: Hệ số góc đối với biến độc lập; + X: Biến độc lập - đường kính lớn (mm); + : Sai số ngẫu nhiên. Xây dựng phương trình hồi quy tuyến tính với biến phụ thuộc Y là khối lượng trứng và biến độc lập X là đường kính lớn được thực hiện như sau: Stat ➔ Regression ➔ Regression... 113
  4. Khai báo các biến khối lượng trứng (P) vào ô Response. Khai báo các biến đường kính lớn (DKL) vào ô Predictors. Chọn OK để có kết quả. Kết quả out put trong cửa sổ Session của Minitab: Regression Analysis: P versus DKL  The regression equation is  P = - 53.7 + 2.04 DKL Predictor Coef SE Coef T P Constant -53.67 12.78 -4.20 0.000 DKL 2.0379 0.2250 9.06 0.000 S = 2.69651 R-Sq = 80.4% R-Sq(adj) = 79.4% Analysis of Variance Source DF SS MS F P Regression 1 596.60 596.60 82.05 0.000 Residual Error 20 145.42 7.27 Total 21 742.02 Unusual Observations Obs DKL P Fit SE Fit Residual St Resid 7 57.1 67.900 62.629 0.579 5.271 2.00R 8 58.2 59.000 64.871 0.658 -5.871 -2.25R R denotes an observation with a large standardized residual. Trong đó:  Kết quả phân tích hồi quy;  Phương trình hồi quy với biến phụ thuộc là khối lượng (Y) và biến độc lập là đường kính lớn (X) từ phần mềm Minitab thu được như sau: Khối lượng trứng = - 53,67124 + 2,03786 × đường kính lớn; 114
  5.  Các tham số ước tính và  Sai số tiêu chuẩn tương ứng đối với từng hệ số của phương trình hồi quy;  Giá trị xác suất của tung độ gốc (P = 0,000) và hệ số góc (P = 0,000) đều < 0,001, nên kết luận các hệ số trong phương trình hồi quy khác 0 (P
  6. Khai báo các biến khối lượng trứng (P) vào ô Response. Khai báo các biến đường kính lớn (DKL) và đường kính nhỏ (DKN) vào ô Predictors. Chọn OK để có kết quả. Kết quả out put trong cửa sổ Session của Minitab: Regression Analysis: P versus DKL, DKN  The regression equation is  P = - 117 + 1.21 DKL + 2.48 DKN Predictor Coef SE Coef T P Constant -116.555 5.472 -21.30 0.000 DKL 1.21473 0.08323 14.60 0.000 DKN 2.4764 0.1623 15.26 0.000 S = 0.759757 R-Sq = 98.5% R-Sq(adj) = 98.4% Analysis of Variance Source DF SS MS F P Regression 2 731.05 365.53 633.24 0.000 Residual Error 19 10.97 0.58 Total 21 742.02 Source DF Seq SS DKL 1 596.60 DKN 1 134.46 Trong đó:  Kết quả phân tích hồi quy;  Phương trình hồi quy với biến phụ thuộc là khối lượng (Y) và biến độc lập là đường kính lớn (X) từ phần mềm Minitab thu được như sau: Khối lượng trứng = -117 + 1,21 × Đường kính lớn + 2,48 × đường kính nhỏ;  Các tham số ước tính và  Sai số tiêu chuẩn tương ứng đối với từng hệ số của phương trình hồi quy;  Giá trị xác suất P = 0,000 ta có thể kết luận các hệ số trong phương trình hồi quy khác 0 (P
  7. trình đóng góp không đáng kể cho mô hình, vì vậy cần loại bỏ khi xây dựng phương trình. Để có thể tìm được phương trình hồi quy tốt nhất cần có thêm các lệnh phụ khác. 3.2.2.1. Phương pháp Best Subset Đây là phương pháp liệt kê tất cả các trường hợp (subset) có thể của các biến độc lập trong mô hình thống kê kèm theo các thông tin (R2 CP và S) để lựa chọn subset tốt nhất. Để đánh giá mức độ đóng góp của biến độc lập X1 là đường kính lớn và X2 là đường kính nhỏ với biến phụ thuộc Y là khối lượng trứng được thực hiện như sau: Stat ➔ Regression ➔ Best Subsets... Khai báo các biến khối lượng trứng (P) vào ô Response. Khai báo các biến đường kính lớn (DKL) và đường kính nhỏ (DKN) vào ô Free predictors. Chọn OK để có kết quả. Kết quả out put trong cửa sổ Session của Minitab: Best Subsets Regression: P versus DKL, DKN  Response is P D D Mallows K K Vars R-Sq R-Sq(adj) Cp S L N 1 82.0 81.0 214.0 2.5878 X 1 80.4 79.4 233.9 2.6965 X 2 98.5 98.4 3.0 0.75976 X X Trong đó:  Kết quả phân tích hồi quy theo phương pháp Best Subset;  Số biến độc lập trong mô hình hồi quy;  Hệ số xác định R² và  Hệ số xác định hiệu chỉnh. Hệ số xác định và hệ số xác định hiệu chỉnh lớn nhất là mô hình tối ưu nhất; 117
  8. + Với 1 biến đường kính nhỏ (DKN) giải thích được 82% mức độ biến thiên của biến khối lượng trứng; + Với 1 biến đường kính lớn (DKL) giải thích được 80,4% mức độ biến thiên của biến khối lượng trứng; + Với 2 biến đường kính lớn (DKL) và đường kính nhỏ (DKN) giải thích được 98,5% mức độ biến thiên của biến khối lượng trứng.  Hệ số Cp, mô hình có hệ số Cp thấp nhất là mô hình tốt nhất nhưng không phải là mô hình tối ưu nhất; CP: Conceptual predictive criterion (còn gọi là hệ số Mallow) được sử dụng để xác định mô hình có khả năng giải thích sự biến động lớn nhất với số biến ít nhất có thể và được tính theo công thức Cp = p + (MSE −  ² )(n − p ) , trong đó ² MSE = trung bình bình phương sai số ngẫu nhiên của mô hình rút gọn, ² = phương sai của mô hình đầy đủ, n = dung lượng mẫu, p = tham số của mô hình rút gọn, mô hình tốt nhất khi có CP nhỏ và gần bằng p;  Giá trị sai số tiêu chuẩn (S) của mô hình hồi quy. Giá trị S nhỏ nhất là mô hình tối ưu nhất;  Tên các biến độc lập trong mô hình. + Với 2 biến đường kính lớn (DKL) và đường kính nhỏ (DKN): Sai số tiêu chuẩn của mô hình bằng 0,76 (nhỏ nhất trong 3 mô hình). Do đó, mô hình gồm hai biến DKL và DKN là tối ưu nhất. 3.2.2.2. Phương pháp Stepwise Để đánh giá mức độ đóng góp của biến độc lập X1 là đường kính lớn và X2 là đường kính nhỏ với biến phụ thuộc Y là khối lượng trứng, ta thực hiện như sau: Stat ➔ Regression ➔ Stepwise... Khai báo các biến khối lượng trứng (P) vào ô Response. Khai báo các biến đường kính lớn (DKL) và đường kính nhỏ (DKN) vào ô Predictors. 118
  9. Chọn Option Methods để lựa chọn phương pháp tính (Forward hoặc Backward). Minitab mặc định sử dụng giá trị alpha (Use alpha values = 0,15, tức biến có giá trị xác suất P < 0,15 mới có thể thêm vào mô hình hồi quy) và phương pháp Forward và Backward (Stepwise). Tuy nhiên, sử dụng giá trị F để thay thế cho giá trị alpha (Minitab mặc định giá trị F = 4) hoặc có thể sử dụng Forward selection hoặc Backward elimination. Chọn OK để có kết quả. Kết quả out put trong cửa sổ Session của Minitab: Stepwise Regression: P versus DKL, DKN  Alpha-to-Enter: 0.15 Alpha-to-Remove: 0.15 Response is P on 2 predictors, with N = 22 Step 1 2 Constant -115.5 -116.6 DKN 4.01 2.48 T-Value 9.53 15.26 P-Value 0.000 0.000 DKL 1.215 T-Value 14.60 P-Value 0.000 S 2.59 0.760 R-Sq 81.95 98.52 R-Sq(adj) 81.05 98.37 Mallows Cp 214.0 3.0 Trong đó:  Kết quả phân tích hồi quy bằng phương pháp Stepwise đối với biến khối lượng (P);  Giá trị xác suất alpha để thêm vào hoặc loại bỏ biến ra khỏi mô hình hồi quy (alpha = 0,15);  Bước 1 (Step 1);  Đưa biến độc lập thứ nhất (DKN) vào trong 119
  10. mô hình với thông tin thu được R² = 81,95% và Cp = 214,0; Các tham số ước tính ở Step 1 (tọa độ gốc -115,5 và hệ số góc đối với biến độc lập DKN = 4,01);  Bước 2 (Step 2);  Đưa thêm biến độc lập thứ hai (DKL) vào trong mô hình với thông tin thu được R² = 98,52% và Cp = 3,0; Các tham số ước tính ở Step 2 (tọa độ gốc -116,6 và hệ số góc đối với biến độc lập DKN = 2,48 và DL = 1,215). - Forward Selection (FORWARD) là phương pháp mà các biến độc lập lần lượt được đưa vào mô hình, mỗi lần đưa một biến vào mô hình gọi là một bước (step). Biến độc lập có giá trị xác suất P
  11. mô hình hồi quy (alpha = 0,25);  Bước 1 (Step 1);  Đưa biến độc lập thứ nhất (DKN) vào trong mô hình với thông tin thu được R² = 81,95% và Cp = 214,0; Các tham số ước tính ở Step 1 (tọa độ gốc -115,5 và hệ số góc đối với biến độc lập DKN = 4,01);  Bước 2 (Step 2);  Đưa thêm biến độc lập thứ hai (DKL) vào trong mô hình với thông tin thu được R² = 98,52% và Cp = 3,0; Các tham số ước tính ở Step 2 (tọa độ gốc -116,6 và hệ số góc đối với biến độc lập DKN = 2,48 và DL = 1,215). - Backward Elimination (BACKWARD) là phương pháp loại bỏ dần những biến độc lập trong phương trình hồi quy (ngược lại với phương pháp Forward Selection). Biến độc lập có giá trị xác suất P >  sẽ bị loại bỏ khỏi mô hình. Chọn Option Methods để lựa chọn phương pháp tính (Backward). Stepwise Regression: P versus DKL, DKN  Backward elimination. Alpha-to-Remove: 0.1 Response is P on 2 predictors, with N = 22 Step 1 Constant -116.6 DKL 1.215 T-Value 14.60 P-Value 0.000 DKN 2.48 T-Value 15.26 P-Value 0.000 S 0.760 R-Sq 98.52 R-Sq(adj) 98.37 Mallows Cp 3.0 121
  12. Trong đó:  Kết quả phân tích hồi quy với phương pháp Stepwise (Backward) đối với biến khối lượng (P);  Giá trị xác suất alpha để thêm vào hoặc loại bỏ biến ra khỏi mô hình hồi quy (alpha = 0,1);  Bước 1 (Step 1);  Đưa tất cả các biến độc lập vào trong mô hình với thông tin thu được R² = 98,52% và Cp = 3,0; Các tham số ước tính ở Step 1 (tọa độ gốc -116,6 và hệ số góc đối với biến độc lập DKL = 1,215 và DKN = 2,48);  Không biến nào bị loại khỏi mô hình, tất cả các biến độc lập (DKL, DKN) đều có trong mô hình cuối cùng ở mức P < 0,15. Phương pháp Stepwise chỉ dựa trên phân tích phần dư trong phân tích hồi quy bội tuyến tính để lựa chọn biến đưa vào và biến đưa ra. Để đảm bảo trong mô hình cuối cùng có mặt một số biến có ý nghĩa về mặt chuyên môn các nhà thống kê thường đưa thêm phần chọn biến bắt buộc có mặt trước khi bắt đầu việc lọc theo Stepwise. 3.2.3. Hồi quy đa thức Mối liên hệ tuyến tính giữa các biến đã được đề cập ở phần trên. Trong sinh học, mối liên hệ giữa các biến có thể tồn tại nhưng ở dưới dạng phi tuyến. Trong những trường hợp như vậy, mô hình tuyến tính nêu trên sẽ không còn phù hợp. Ví dụ 3.2: Khối lượng của gà Hồ (g) qua các tuần tuổi (từ 0 đến 10 tuần tuổi) được trình bày ở bảng sau: Tuần tuổi (tuần) 0 1 2 3 4 5 6 7 8 9 10 Khối lượng (g) 35 45 72 134 199 303 374 529 672 717 899 Nếu thể hiện mối liên hệ giữa khối lượng và tuần tuổi bằng các điểm () trên đồ thị, mối liên hệ giữa 2 biến này không phải là tuyến tính (đường ---) mà có mối liên hệ phi tuyến bậc 2 (đường —). 122
  13. 1000 y = 6.6119x2 + 9.5028x + 0.5636 R² = 0.9931 800 600 Khối lượng (g) y = 88.845x - 171.35 R² = 0.952 400 200 0 0 1 2 3 4 5 6 7 8 9 10 Tuần tuổi -200 Mô hình phân tích có dạng: y = a + b1x + b2x² Cấu trúc số liệu Số liệu được nhập vào trong cửa sổ Worksheet của phần mềm Minitab; bao gồm 2 cột: 1) Cột Tuần tuổi C1 (TT); 2) Cột Khối lượng C2 (KL) . Tạo cột TT^2 C5 (TT2) Calc ➔ Calculator... 123
  14. Khai báo biến TT2 vào ô Store result in variable. Nhập công thức: TT*TT vào ô Expression. Chọn OK để có kết quả. Stat ➔ Regression ➔ Regression... 124
  15. Khai báo biến KL vào ô Response. Khai báo biến TT và TT2 vào ô Predictors. Chọn OK để có kết quả. Regression Analysis: KL versus TT, TT2  The regression equation is KL = 16.7 + 22.7 TT + 6.61 TT2 Predictor Coef SE Coef T P Constant 16.68 21.35 0.78 0.457 TT 22.727 9.933 2.29 0.051 TT2 6.6119 0.9567 6.91 0.000 S = 28.0230 R-Sq = 99.3% R-Sq(adj) = 99.1% Analysis of Variance Source DF SS MS F P Regression 2 905796 452898 576.73 0.000 Residual Error 8 6282 785 Total 10 912078 Source DF Seq SS TT 1 868287 TT2 1 37509 trong đó:  Kết quả phân tích hồi quy đối với biến phụ thuộc khối lượng (KL);  Mô hình hồi quy có thể được viết như sau: Khối lượng = 16,68 + 22,7 × Tuần + 6,61 × Tuần²  Tung độ gốc (Intercept) =16,68;  Hệ số góc đối với biến thứ nhất (TT) = 22,7 và  Hệ số góc đối với biến thứ 2 (TT2) = 6,61;  Xác suất đối với từng hệ số;  Hệ số xác định của mô hình (R-Sq) = 99,3% (nghĩa là 99,31% biến động được giải thích bằng mô 125
  16. hình hồi quy);  Hệ số xác định hiệu chỉnh theo các biến độc lập trong mô hình (Adj R- Sq) = 99,1%.  Adj R-Sq = 1 − (1 − R ) ( N − 1)  , trong đó R², N và p lần lượt là hệ số xác định, dung 2   N − p −1    lượng mẫu và số biến độc lập;  Bảng kết quả phân tích phương sai. 3.2.4. Hồi quy phi tuyến Ví dụ 3.3: Khối lượng của gà Hồ (g) qua các tuần tuổi (từ 0 đến 32 tuần tuổi) được trình bày ở bảng sau: Tuần tuổi 0 1 2 3 4 5 6 7 8 9 10 Khối lượng 35 45 72 134 199 303 374 529 672 717 899 Tuần tuổi 11 12 13 14 15 16 17 18 19 20 21 Khối lượng 933 1042 1228 1386 1541 1729 1711 1726 1771 1927 2032 Tuần tuổi 22 23 24 25 26 27 28 29 30 31 32 Khối lượng 2158 2232 2557 2632 2718 2818 2937 2960 2964 2967 2977 126
  17. Số liệu có thể phân tích theo hàm sinh trưởng Gompertz theo mô hình sau: −𝒃𝒙 Y = 𝒎𝒆−𝒂𝒆 = m * EXP [-a * EXP (-bx)] Trong đó: + Y = khối lượng gà (g); + m = khối lượng tiệm cận trên – khối lượng trưởng thành (g); + a = hằng số tích hợp liên quan đến khối lượng sơ sinh (0 tuần tuổi); + b = tỷ lệ tốc độ tăng trưởng tối đa so với khối lượng trưởng thành của gà; + x = tuổi (tuần); + EXP = cơ số logarit tự nhiên e (2,71828). Cấu trúc số liệu 127
  18. Số liệu được nhập vào trong cửa sổ Worksheet của phần mềm Minitab; bao gồm 2 cột: 1) Cột Tuần tuổi C1 (TT); 2) Cột Khối lượng C2 (KL). Stat ➔ Regression ➔ Nonlinear Regression... Khai báo biến KL vào ô Response. 128
  19. Chọn option Use Calculator của mục Expectation Funtion để khai báo mô hình hồi quy (ví dụ: Gompertz Growth). Chọn option Parameters để khai các tham số tham số Theta1, Theta2 và Theta3. Chọn OK để có kết quả. Nonlinear Regression: KL = Theta1 * exp(-Theta2 * exp(-Theta3 * TT))  129
  20. Method Algorithm Gauss-Newton Max iterations 200 Tolerance 0.00001 Starting Values for Parameters Parameter Value Theta1 2980 Theta2 4.4 Theta3 0.14 Equation KL = 3734.33 * exp(-3.86255 * exp(-0.0932028 * TT)) Parameter Estimates Parameter Estimate SE Estimate Theta1 3734.33 155.331 Theta2 3.86 0.191 Theta3 0.09 0.006 KL = Theta1 * exp(-Theta2 * exp(-Theta3 * TT)) Lack of Fit There are no replicates. Minitab cannot do the lack of fit test based on pure error. Summary Iterations 8 Final SSE 239613 DFE 30 MSE 7987.09 S 89.3705 Trong đó:  Kết phân tích hồi quy phi tuyến đối với biến phụ thuộc khối lượng (KL);  Giá trị xuất phát «Starting Value», các giá trị này đưa vào mô hình dựa trên số liệu của mô hình và/hoặc các nghiên cứu tương tự. Ví dụ giá trị m = 2.980 dựa vào khối lượng của gà Hồ ở tuần tuổi 32 (2.977g), giá trị a và b lấy từ kết quả nghiên cứu của (Nassim & cs., 2011);  Mô hình hồi quy phi tuyến để ước tính khối lượng gà Hồ như sau: Khối lượng = 3734,33 × exp(-3,86 × exp(-0,09 × tuần))  Các tham số ước tính (Theta1 = 3734,33; Theta2 = 3,86 và Theta3 = 0,09) và sai số tiêu chuẩn ước tính (SE Estimate); Tuổi tại điểm uốn xPI = ln(a)/b = ln(3,8625)/0,0932 = 14,5 tuần Khối lượng gà Hồ tại điểm uốn YPI = m/e = 3734;3/2,71828 = 1373,77g Tăng khối lượng tuyệt đối (g/tuần) tại điểm uốn: MWGPI = mb/e = 3734.3 × 0.0932/2,71828 = 128g. 130
ADSENSE

CÓ THỂ BẠN MUỐN DOWNLOAD

 

Đồng bộ tài khoản
2=>2