Bình phương nhỏ nhất có trọng số python scipy

Có vẻ như vấn đề là toán học, chứ không phải mã hóa. Tôi sẽ không giúp với con trăn, vì bạn biết rõ mình đang làm gì, nhưng tôi chắc rằng mọi người sẽ đánh giá cao việc đăng lại giải pháp của bạn

Tôi có thể nghĩ ra hai phương pháp, một phương pháp thô sơ và nhanh chóng, phương pháp kia nghiêm ngặt hơn

Phương pháp 1.
- Tạo trọng số số nguyên, nhưng đảo ngược sai số [1/lỗi], nhân với một hằng số phù hợp và làm tròn thành số nguyên gần nhất.
- Tạo tập dữ liệu mới bằng cách thêm nhiều bản sao của mỗi điểm dữ liệu, tương ứng với số nguyên trên.
- Thực hiện điều chỉnh tối thiểu bình phương trên tập dữ liệu mới này.

Rõ ràng bằng cách chọn hằng số lớn phù hợp, bạn có thể lấy trọng số khá chính xác. Ưu điểm lớn là đó là một tinh chỉnh nhỏ trên mã của bạn

Đây là sự tái tạo của khóa học Stanford Stats 191 [xem https. // trang web. standford. edu/class/stats191/], sử dụng các công cụ hệ sinh thái Python, thay vì R. Đây là bài giảng "Phép đổi và Bình phương nhỏ nhất có trọng số"

Thiết lập Sổ tay Ban đầu¶

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
6 ghi lại môi trường Python và gói hiện tại,
fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
7 là trình định dạng Python ưa thích của tôi

Trong 2]

%load_ext watermark

Trong 3]

%load_ext lab_black

Trong [4]

%matplotlib inline

Tất cả quá trình nhập đều có ở đây [không phải tất cả đều được sử dụng trong Notebook này]

Trong [43]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings

Mô hình phi tuyến¶

Chúng tôi xem xét một tập dữ liệu trong đó mô hình tuyến tính không phù hợp

Tải và Khám phá Tập dữ liệu¶

Trong [6]

data = pd.read_csv['../data/bacteria.txt', sep='\t']

Trong [7]

data.head[]

Ra[7]

tN_t0135512211231973416645142

Trong [8]

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']

Ra[8]

Text[0, 0.5, 'Count']

Điều chỉnh mô hình tuyến tính¶

Trong [9]

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
0

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
1

Ra[9]

Kết quả hồi quy OLSDep. Biến đổi. N_tR bình phương. 0. 823Mô hình. OLSAdj. bình phương R. 0. 810Phương pháp. Số bình phương nhỏ nhất F-statistic. 60. 62Ngày. Thứ ba, ngày 07 tháng 4 năm 2020Prob [F-statistic]. 3. 01e-06Thời gian. 16. 42. 55Log-Khả năng. -76. 216Không. quan sát. 15AIC. 156. Dư 4Df. 13BICK. 157. Mẫu 8Df. 1Loại hiệp phương sai. nonrobustcoefstd errtP>. t. [0. 0250. 975]Đánh chặn259. 581022. 73011. 4200. 000210. 476308. 686t-19. 46432. 500-7. 7860. 000-24. 865-14. 063toàn bộ. 14. 123Durbin-Watson. 0. 803Prob[Xe đa năng]. 0. 001Jarque-Bera [JB]. 10. 490Xiên. 1. 648Prob[JB]. 0. 00527Kurtosis. 5. 432Cond. Không. 19. 3

Cảnh báo.
[1] Lỗi chuẩn giả định rằng ma trận hiệp phương sai của các lỗi được chỉ định chính xác.

Trực quan hóa sự phù hợp của mô hình tuyến tính của chúng tôi với dữ liệu

Trong [10]

%load_ext lab_black
0

Ra[10]

%load_ext lab_black
1

Chẩn đoán hồi quy¶

Các ghi chú bài giảng STATS191 có bốn biểu đồ có thể được sử dụng để đánh giá chất lượng mô hình tuyến tính. Chúng tôi xác định một hàm để hiển thị những điều này, được cung cấp đối tượng RegressionResults từ một mô hình phù hợp với

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
8

Trong [77]

%load_ext lab_black
2

Trong [76]

%load_ext lab_black
3

Chúng tôi thấy rằng phần dư rõ ràng không được phân phối đồng đều [Residuals vs Fitted graph]

Hơn nữa, nếu chúng ta nhìn vào đồ thị ảnh hưởng và đòn bẩy của

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
9, chúng ta sẽ thấy một ngoại lệ rõ ràng

Trong [13]

%load_ext lab_black
4

Trong [14]

%load_ext lab_black
5

Phù hợp với số mũ¶

Ở giai đoạn này, chúng tôi quyết định rằng một mô hình tốt hơn có thể là mô hình phân rã theo hàm mũ, trong đó nhật ký của

Text[0, 0.5, 'Count']
0 có thể là tuyến tính theo thời gian

Lưu ý rằng chúng ta có thể sử dụng các hàm

Text[0, 0.5, 'Count']
1 bên trong công thức
fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
8

Trong [15]

%load_ext lab_black
6

fig, ax = plt.subplots[figsize=[8, 8]]

ax.plot[data['t'], data['N_t'], 'd']
ax.set_xlabel['Time']
ax.set_ylabel['Count']
1

Ra[15]

Kết quả hồi quy OLSDep. Biến đổi. np. log[N_t]R bình phương. 0. 988Mô hình. OLSAdj. bình phương R. 0. 987Phương pháp. Số bình phương nhỏ nhất F-statistic. 1104. Ngày tháng. Thứ ba, ngày 07 tháng 4 năm 2020Prob [F-statistic]. 5. 86e-14Thời gian. 16. 42. 57Log-Khả năng. 12. 896Không. quan sát. 15AIC. -21. Dư lượng 79Df. 13BICK. -20. Mẫu 38Df. 1Loại hiệp phương sai. nonrobustcoefstd errtP>. t. [0. 0250. 975]Đánh chặn5. 97320. 06099. 9220. 0005. 8446. 102t-0. 21840. 007-33. 2220. 000-0. 233-0. 204Xe vạn năng. 0. 024Durbin-Watson. 2. 659Prob[Xe đa năng]. 0. 988Jarque-Bera [JB]. 0. 204Xiên. -0. 074Prob[JB]. 0. 903Kiểu nhọn. 2. 448Điều kiện. Không. 19. 3

Cảnh báo.
[1] Lỗi chuẩn giả định rằng ma trận hiệp phương sai của các lỗi được chỉ định chính xác.

Đánh giá mô hình hàm mũ¶

Vẽ sơ đồ các dự đoán của mô hình so với thực tế theo hai cách khác nhau, chúng ta có thể thấy chúng ta có sự phù hợp tốt hơn nhiều [đồng thời giá trị R^2 gần với 1 hơn nhiều. 0]

Trong [16]

%load_ext lab_black
8

Ra[16]

Text[0, 0.5, 'Count']

Trong [17]

%matplotlib inline
0

%matplotlib inline
1

Ra[17]

%matplotlib inline
2

Đồ thị chẩn đoán cũng trông đẹp hơn [đặc biệt là Phần dư so với Đã trang bị]

Trong [18]

%matplotlib inline
3

Trong 19]

%matplotlib inline
4

Trong 20]

%matplotlib inline
5

Phương sai không cố định¶

Bây giờ chúng ta xem xét một tập dữ liệu, với một thuật ngữ lỗi không cố định từ quan sát này sang quan sát khác [xem bên dưới để biết mô tả về tập dữ liệu]

Mô tả biến

Y. Chi tiêu giáo dục bình quân đầu người theo nhà nước X1. Thu nhập bình quân đầu người năm 1973 theo bang X2. Tỷ trọng dân số dưới 18 tuổi X3. Tỷ trọng khu vực thành thị Vùng. Các bang nằm ở khu vực nào của đất nước

Đọc và khám phá tập dữ liệu¶

Trong [21]

%matplotlib inline
6

Trong [22]

%matplotlib inline
7

Hết[22]

STATEYX1X2X3Region0ME235394432550811NH231457832356412VT270401132832213MA261523330584614RI30047803038711

Thật khó chịu, chúng tôi có dấu cách trong tên cột, vì vậy hãy sửa chúng

Trong [23]

%matplotlib inline
8

Hết[23]

%matplotlib inline
9

Trong [24]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
0

Ra[24]

STATEYX1X2X3Region0ME23539443255081

Trong [25]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
1

Hết[25]

YX1X2X3Regioncount50. 00000050. 00000050. 00000050. 00000050. 000000mean284. 6000004675. 120000325. 740000657. 8000002. 660000std61. 340136644. 50625419. 423119145. 0163961. 061574min208. 0000003448. 000000287. 000000322. 0000001. 00000025%234. 2500004137. 250000310. 750000546. 7500002. 00000050%269. 5000004706. 000000324. 500000662. 5000003. 00000075%316. 7500005054. 250000333. 000000782. 2500003. 750000max546. 0000005889. 000000386. 000000909. 0000004. 000000

Fit Linear Model¶

Trong [26]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
2

Hết[26]

Kết quả hồi quy OLSDep. Biến đổi. YR bình phương. 0. 591Mô hình. OLSAdj. bình phương R. 0. Phương pháp 565. Số bình phương nhỏ nhất F-statistic. 22. 19Ngày. Thứ ba, ngày 07 tháng 4 năm 2020Prob [F-statistic]. 4. 94e-09Thời gian. 16. 43. 02Khả năng đăng nhập. -253. 89Không. quan sát. 50AIC. 515. Dư 8Df. 46BIC. 523. Mô hình 4Df. 3Loại hiệp phương sai. nonrobustcoefstd errtP>. t. [0. 0250. 975]Đánh chặn-556. 5680123. 195-4. 5180. 000-804. 547-308. 589X10. 07240. 0126. 2390. 0000. 0490. 096X21. 55210. 3154. 9320. 0000. 9192. 185X3-0. 00430. 051-0. 0830. 934-0. 1080. 099Toàn diện. 0. 719Durbin-Watson. 2. 149Prob[Xe đa năng]. 0. 698Jarque-Bera [JB]. 0. 760Xiên. 0. 264Prob[JB]. 0. 684Kurtosis. 2. 705Cond. Không. 1. 03e+05

Cảnh báo.
[1] Lỗi chuẩn giả định rằng ma trận hiệp phương sai của các lỗi được chỉ định chính xác.
[2] Số điều kiện lớn, 1. 03e+05. Điều này có thể chỉ ra rằng có
đa cộng tuyến mạnh hoặc các vấn đề về số khác.

Đánh giá mô hình¶

Biểu đồ Dư so với Đã trang bị cho thấy mức tăng rõ ràng về giá trị còn lại khi giá trị được trang bị tăng [i. e. phương sai của thuật ngữ lỗi không phải là hằng số]

Trong [27]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
3

Hơn nữa, chúng tôi dường như có một quan sát ngoại lệ [Bang Alaska, mã AK]

Trong [28]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
4

Trong [29]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
5

Hết[29]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
6

Trong [30]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
7

Bây giờ chúng ta trực quan hóa phần dư [được chuẩn hóa] theo vùng và có vẻ như mỗi vùng có một phương sai khác nhau của thuật ngữ lỗi

Trong [78]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
8

Hết[78]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
9

Thật khó chịu, chúng tôi thấy các mã trạng thái cũng có dấu cách, vì vậy hãy sửa các mã này

Trong [34]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
0

Ra[34]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
1

Trong [35]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
2

Ra[35]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
3

Bây giờ chúng tôi làm sạch tập dữ liệu, bằng cách loại trừ tiểu bang Alaska và tạo một bản sao của DataFrame

Sau đó chúng tôi khớp lại một mô hình tuyến tính

Trong [36]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
4

data = pd.read_csv['../data/bacteria.txt', sep='\t']
5

Ra[36]

Kết quả hồi quy OLSDep. Biến đổi. YR bình phương. 0. 497Mô hình. OLSAdj. bình phương R. 0. 463phương pháp. Số bình phương nhỏ nhất F-statistic. 14. 80Ngày. Thứ ba, ngày 07 tháng 4 năm 2020Prob [F-statistic]. 7. 65e-07Thời gian. 17. 05. 11Khả năng đăng nhập. -242. 77Không. quan sát. 49AIC. 493. Dư 5Df. 45BIC. 501. Mô hình 1Df. 3Loại hiệp phương sai. nonrobustcoefstd errtP>. t. [0. 0250. 975]Đánh chặn-277. 5773132. 423-2. 0960. 042-544. 291-10. 864X10. 04830. 0123. 9760. 0000. 0240. 073X20. 88690. 3312. 6780. 0100. 2201. 554X30. 06680. 0491. 3540. 183-0. 0330. 166Xe buýt. 0. 883Durbin-Watson. 1. 946Prob[Xe buýt đa năng]. 0. 643Jarque-Bera [JB]. 0. 932Xiên. 0. 289Prob[JB]. 0. 628Kurtosis. 2. 652Cond. Không. 1. 23e+05

Cảnh báo.
[1] Lỗi chuẩn giả định rằng ma trận hiệp phương sai của các lỗi được chỉ định chính xác.
[2] Số điều kiện lớn, 1. 23e+05. Điều này có thể chỉ ra rằng có
đa cộng tuyến mạnh hoặc các vấn đề về số khác.

Đánh giá mô hình¶

Làm sạch dữ liệu không khắc phục được phần dư phương sai non_constant của chúng tôi

Trong [39]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
6

Dư địa theo Vùng vẫn chênh lệch giữa các Vùng

Trong [40]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
7

Ra[40]

import pandas as pd
import numpy as np
import seaborn as sns

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

from scipy import stats
from statsmodels.formula.api import ols
from statsmodels.formula.api import rlm
from statsmodels.formula.api import mixedlm
from statsmodels.formula.api import wls
import statsmodels.api as sm

from patsy import dmatrices

from patsy.contrasts import Treatment
from patsy.contrasts import Sum
from patsy.contrasts import Diff


import warnings
9

Bình phương nhỏ nhất có trọng số¶

Chúng tôi xác định một mô hình mới, trong đó chúng tôi đã chuẩn hóa phương sai bằng trọng số với giá trị của Y [giả sử một mô hình có số hạng sai số là N[0,eps * Y]. Chúng tôi giả sử các quan sát trong mỗi Vùng có số hạng sai số phương sai không đổi

Mảng

Text[0, 0.5, 'Count']
3 sẽ giữ [1/phương sai] của thuật ngữ lỗi [trong mô hình của chúng tôi]

Trong [47]

data = pd.read_csv['../data/bacteria.txt', sep='\t']
9

Nhận danh sách mã Vùng duy nhất

Trong [50]

data.head[]
0

Hết[50]

data.head[]
1

Đối với mỗi Khu vực, hãy tính ước tính phương sai của thuật ngữ sai số cho khu vực đó và điền trọng số cho tất cả các quan sát từ khu vực đó

Trong [52]

data.head[]
2

data.head[]
3

Lưu ý rằng giá trị của các trọng số ở trên phù hợp với ghi chú bài giảng STATS191

Làm vừa người mẫu với quả nặng¶

Trong [59]

data.head[]
4

Hết[59]

Kết quả hồi quy WLSDep. Biến đổi. YR bình phương. 0. 757Mô hình. WLSAdj. bình phương R. 0. 741Phương pháp. Số bình phương nhỏ nhất F-statistic. 46. 75Ngày. Thứ ba, ngày 07 tháng 4 năm 2020Prob [F-statistic]. 7. 10e-14Thời gian. 18. 22. 54Log-Khả năng. -231. 58Không. quan sát. 49AIC. 471. Dư 2Df. 45BIC. 478. Mẫu 7Df. 3Loại hiệp phương sai. nonrobustcoefstd errtP>. t. [0. 0250. 975]Đánh chặn-315. 517078. 178-4. 0360. 000-472. 975-158. 059X10. 06230. 0087. 9160. 0000. 0460. 078X20. 87430. 2004. 3660. 0000. 4711. 278X30. 02940. 0340. 8570. 396-0. 0400. 098Xe vạn năng. 6. 395Durbin-Watson. 2. 023Prob[Xe đa năng]. 0. 041Jarque-Bera [JB]. 2. 405Xiên. 0. 147Prob[JB]. 0. 300Kurtosis. 1. 955Điều kiện. Không. 1. 05e+05

Cảnh báo.
[1] Lỗi chuẩn giả định rằng ma trận hiệp phương sai của các lỗi được chỉ định chính xác.
[2] Số điều kiện lớn, 1. 05e+05. Điều này có thể chỉ ra rằng có
đa cộng tuyến mạnh hoặc các vấn đề về số khác.

Các kết quả trên phù hợp với kết quả STATS191 ở hai chữ số thập phân đầu tiên [tôi nghi ngờ cảnh báo

Số điều kiện lớn, 1. 05e+05. Điều này có thể chỉ ra rằng có hiện tượng đa cộng tuyến mạnh hoặc các vấn đề về số khác

nên được chú ý

Để kiểm tra trọng số, bạn có thể lấy thông qua thuộc tính

Text[0, 0.5, 'Count']
4 của RegressionResults

Trong [79]

data.head[]
5

Hết[79]

data.head[]
6

Có vẻ như sau khi điều chỉnh mô hình Bình phương nhỏ nhất có trọng số, Phần dư Pearson là cách tốt hơn để chuẩn hóa phần dư. Biểu đồ hình hộp bên dưới cho thấy chúng tôi đã tính đến phần lớn sự khác biệt giữa các Khu vực trong mô hình của mình, nhưng không phải tất cả

Chủ Đề