|
688 | 688 | "id": "d0c00cea-1826-43e8-bb1b-45b03af7bf4c",
|
689 | 689 | "metadata": {},
|
690 | 690 | "source": [
|
691 |
| - "## Bonus 1: adding noise" |
| 691 | + "## Adding noise" |
692 | 692 | ]
|
693 | 693 | },
|
694 | 694 | {
|
|
820 | 820 | "The EOF analysis has again successfully identified our two modes. The third and fourth modes explain more variance than in the clean example, but the spatial patterns strongly suggest that they are just noise."
|
821 | 821 | ]
|
822 | 822 | },
|
823 |
| - { |
824 |
| - "cell_type": "markdown", |
825 |
| - "id": "3a1d8fb0-5ef6-461c-a3ab-30061755ba53", |
826 |
| - "metadata": {}, |
827 |
| - "source": [ |
828 |
| - "## Bonus 2: method of snapshots" |
829 |
| - ] |
830 |
| - }, |
831 |
| - { |
832 |
| - "cell_type": "markdown", |
833 |
| - "id": "cda271e3-2730-4646-8068-c8ed5a2a1766", |
834 |
| - "metadata": {}, |
835 |
| - "source": [ |
836 |
| - "Sometimes you have a lot more locations than time steps in your data, making computing the covariance matrix impractical, if not impossible. For example, consider a few years of 0.1° gridded climate data formatted as 32-bit floats. That's $(3600\\times1800)^2\\times 4\\approx 1.68\\times10^{14}$ bytes, or about 168 TB, required to hold the covariance matrix in memory.\n", |
837 |
| - "\n", |
838 |
| - "To get around this issue, we turn to the _method of snapshots_ {cite:p}`sirovich_turbulence_1987`. The only change in our simple example is that we compute the covariance of the transposed data, resulting in a matrix of dimension (time, time), which would be much more manageable in our 0.1° data example. Here is the result for our clean 2-mode data:" |
839 |
| - ] |
840 |
| - }, |
841 |
| - { |
842 |
| - "cell_type": "code", |
843 |
| - "execution_count": null, |
844 |
| - "id": "86240329-089d-4fea-8964-efc68428c330", |
845 |
| - "metadata": {}, |
846 |
| - "outputs": [], |
847 |
| - "source": [ |
848 |
| - "s_R = np.cov(data_2d.T)\n", |
849 |
| - "(s_eigval, s_eigvec) = np.linalg.eig(s_R)\n", |
850 |
| - "\n", |
851 |
| - "s_i = np.argsort(-s_eigval)\n", |
852 |
| - "s_eigval = -np.sort(-s_eigval)\n", |
853 |
| - "s_eigvec = s_eigvec[:,s_i]\n", |
854 |
| - "s_eigval = s_eigval[:r]\n", |
855 |
| - "s_eigvec = s_eigvec[:, :r]\n", |
856 |
| - "\n", |
857 |
| - "s_pvar = np.real(s_eigval/np.trace(s_R)*100)\n", |
858 |
| - "s_pcs = np.real(s_eigvec)\n", |
859 |
| - "s_eofs = np.real(np.dot(data_2d, s_eigvec)).reshape(10, 20, r)" |
860 |
| - ] |
861 |
| - }, |
862 |
| - { |
863 |
| - "cell_type": "code", |
864 |
| - "execution_count": null, |
865 |
| - "id": "c9b85632-3170-4803-ac5e-3d7c0f5083bb", |
866 |
| - "metadata": {}, |
867 |
| - "outputs": [], |
868 |
| - "source": [ |
869 |
| - "plt.plot(s_pcs[:, 0], label='PC1')\n", |
870 |
| - "plt.plot(s_pcs[:, 1], label='PC2')\n", |
871 |
| - "plt.plot(s_pcs[:, 2], label='PC3', linestyle=':')\n", |
872 |
| - "plt.plot(s_pcs[:, 3], label='PC4', linestyle=':')\n", |
873 |
| - "plt.legend(loc='upper left')\n", |
874 |
| - "plt.ylabel('PC')\n", |
875 |
| - "plt.xlabel('$t$')\n", |
876 |
| - "plt.xlim(0, nt-1)\n", |
877 |
| - "plt.title('First four PCs (method of snapshots)')" |
878 |
| - ] |
879 |
| - }, |
880 |
| - { |
881 |
| - "cell_type": "code", |
882 |
| - "execution_count": null, |
883 |
| - "id": "e4501a34-a714-494e-a5c0-6ad583a19505", |
884 |
| - "metadata": {}, |
885 |
| - "outputs": [], |
886 |
| - "source": [ |
887 |
| - "s_fig, s_ax = plt.subplots(2, 2, figsize=(10, 5), layout='constrained')\n", |
888 |
| - "\n", |
889 |
| - "mos_eof1_plot = s_ax[0, 0].pcolormesh(x[:, :, 0], y[:, :, 0], s_eofs[:, :, 0], cmap='RdBu_r', norm=CenteredNorm())\n", |
890 |
| - "plt.colorbar(mos_eof1_plot, ax=s_ax[0, 0])\n", |
891 |
| - "s_ax[0, 0].set_ylabel('$y$')\n", |
892 |
| - "s_ax[0, 0].set_title('EOF1 (%.1f%%)' %s_pvar[0])\n", |
893 |
| - "\n", |
894 |
| - "mos_eof2_plot = s_ax[0, 1].pcolormesh(x[:, :, 0], y[:, :, 0], s_eofs[:, :, 1], cmap='RdBu_r', norm=CenteredNorm())\n", |
895 |
| - "plt.colorbar(mos_eof1_plot, ax=s_ax[0, 1])\n", |
896 |
| - "s_ax[0, 1].set_title('EOF2 (%.1f%%)' %s_pvar[1])\n", |
897 |
| - "\n", |
898 |
| - "mos_eof3_plot = s_ax[1, 0].pcolormesh(x[:, :, 0], y[:, :, 0], s_eofs[:, :, 2], cmap='RdBu_r', norm=CenteredNorm())\n", |
899 |
| - "plt.colorbar(mos_eof3_plot, ax=s_ax[1, 0])\n", |
900 |
| - "s_ax[1, 0].set_ylabel('$y$')\n", |
901 |
| - "s_ax[1, 0].set_xlabel('$x$')\n", |
902 |
| - "s_ax[1, 0].set_title('EOF3 (%.1f%%)' %s_pvar[2])\n", |
903 |
| - "\n", |
904 |
| - "mos_eof4_plot = s_ax[1, 1].pcolormesh(x[:, :, 0], y[:, :, 0], s_eofs[:, :, 3], cmap='RdBu_r', norm=CenteredNorm())\n", |
905 |
| - "plt.colorbar(mos_eof4_plot, ax=s_ax[1, 1])\n", |
906 |
| - "s_ax[1, 1].set_xlabel('$x$')\n", |
907 |
| - "s_ax[1, 1].set_title('EOF4 (%.1f%%)' %s_pvar[3])" |
908 |
| - ] |
909 |
| - }, |
910 | 823 | {
|
911 | 824 | "cell_type": "markdown",
|
912 | 825 | "id": "338f68ea-ce77-44df-9a3d-da4d3142e4b3",
|
|
0 commit comments