not much

.gitignore

@@ -184,4 +184,7 @@ ipython_config.py
 
 
 
-.conda/
+.conda*/
+
+data/
+# for rn

.vscode/extensions.json

@@ -0,0 +1,6 @@
+{
+    "recommendations": [
+        "ms-python.python",
+        "ms-toolsai.jupyter"
+    ]
+}

README.md

@@ -14,8 +14,10 @@ how to git? do not `push --force`. prefer pull requests. its like github
 
 we will plan to do this in multiple different ways using scikit-learn and pytorch. tensorflow is ancient.
 
-> [!IMPORTANT]
-> phase out tensorflow. dont use it.
+> [!NOTE]
+> ~phase out tensorflow. dont use it.~
+>
+> its 2 different platforms for 2 different types of people. use whatever ig
 
 maybe also fit more weather related stats? this dataset is taken from nasa so sources are there
 

dataset_update.py

@@ -1,12 +1,45 @@
+import sys
 import urllib.parse
 import urllib.request
 from datetime import datetime
+import logging
+from global_config import DATA_FILE_NAME
+
+# Configure logging
+logging.basicConfig(
+    level=logging.INFO, format="%(asctime)s - %(levelname)s - %(message)s"
+)
 
 # Define parameters
-latitude, longitude = 21.04815477090267, 105.80123159165096  # lmao
-start_date = "2010-01-01"
+latitude, longitude = 21.04815477090267, 105.80123159165096  # Default: Hanoi
+# start_date = "1981-01-01"
+start_date = "2001-01-01"
 end_date = datetime.now().strftime("%Y-%m-%d")  # Today's date
 
+# most used parameters from NASA POWER API:
+# T2M             MERRA-2 Temperature at 2 Meters (C)
+# T2M_MAX         MERRA-2 Maximum Temperature at 2 Meters (C)
+# T2M_MIN         MERRA-2 Minimum Temperature at 2 Meters (C)
+# T2MDEW          MERRA-2 Dew/Frost Point at 2 Meters (C)
+# T2MWET          MERRA-2 Wet Bulb Temperature at 2 Meters (C)
+# PRECTOTCORR     MERRA-2 Precipitation Corrected (mm/day)
+# ALLSKY_SFC_SW_DWN All Sky Surface Shortwave Downward Irradiance (kW-hr/m^2/day)
+# CLRSKY_SFC_SW_DWN Clear Sky Surface Shortwave Downward Irradiance (kW-hr/m^2/day)
+# RH2M            MERRA-2 Relative Humidity at 2 Meters (%)
+# QV2M            MERRA-2 Specific Humidity at 2 Meters (g/kg)
+# PS              MERRA-2 Surface Pressure (kPa)
+# WS10M           MERRA-2 Wind Speed at 10 Meters (m/s)
+# WS50M           MERRA-2 Wind Speed at 50 Meters (m/s)
+# WD10M           MERRA-2 Wind Direction at 10 Meters (Degrees)
+# WD50M           MERRA-2 Wind Direction at 50 Meters (Degrees)
+# TS              MERRA-2 Earth Skin Temperature (C)
+# PRECTOT         MERRA-2 Precipitation (mm/day)
+# ALLSKY_KT       All Sky Insolation Clearness Index (dimensionless)
+# CLRSKY_KT       Clear Sky Insolation Clearness Index (dimensionless)
+# ALLSKY_SFC_LW_DWN All Sky Surface Longwave Downward Irradiance (W/m^2)
+# CLOUD_AMT       Cloud amount (%)
+
+
 # Define desired parameters
 # T2M             MERRA-2 Temperature at 2 Meters (C)
 # PRECTOTCORR     MERRA-2 Precipitation Corrected (mm/day)
@@ -16,9 +49,29 @@ end_date = datetime.now().strftime("%Y-%m-%d")  # Today's date
 # WS10M           MERRA-2 Wind Speed at 10 Meters (m/s)
 # WD10M           MERRA-2 Wind Direction at 10 Meters (Degrees)
 # more at https://power.larc.nasa.gov/data-access-viewer/
+weather_params = [
+    "T2M",
+    "T2M_MAX",
+    "T2M_MIN",
+    "T2MDEW",
+    "PRECTOTCORR",
+    "QV2M",
+    "RH2M",
+    "PS",
+    "CLOUD_AMT",
+    "WS10M",
+    "WD10M",
+    "ALLSKY_KT",
+    "ALLSKY_SFC_SW_DWN",
+]
+
+
+if len(weather_params) > 20:
+    logging.error("too much! %d params bigger than 20", len(weather_params))
+    sys.exit(1)
 
 params = {
-    "parameters": "QV2M,T2M,PRECTOTCORR",
+    "parameters": ",".join(weather_params),
     "community": "RE",
     "longitude": longitude,
     "latitude": latitude,
@@ -37,13 +90,22 @@ url = f"{base_url}?{query_string}"
 
 # Make the API request using urllib
 try:
+    logging.info("Requesting data from NASA POWER API...")
     with urllib.request.urlopen(url) as response:
         if response.status == 200:
-            csv_data = response.read().decode("utf-8")
-            with open("hanoi_weather_data.csv", "w", encoding="utf-8") as file:
-                file.write(csv_data)
-            print("Data successfully downloaded and saved to 'hanoi_weather_data.csv'.")
+            # Save data to file
+            with open(DATA_FILE_NAME, "wb") as file:
+                while chunk := response.read(1024):  # Read 1 KB at a time
+                    file.write(chunk)
+            logging.info(
+                "Data successfully downloaded and saved to '%s'.", DATA_FILE_NAME
+            )
         else:
-            print(f"Failed to retrieve data: HTTP Status Code {response.status}")
+            print("re", response.read().decode())
+            logging.error(
+                "Failed to retrieve data: HTTP Status Code %s", response.status
+            )
 except urllib.error.URLError as e:
-    print(f"Failed to retrieve data: {e.reason}")
+    logging.error("Failed to retrieve data: %s", e)
+except Exception as e:
+    logging.error("An unexpected error occurred: %s", e)

global_config.py

@@ -0,0 +1 @@
+DATA_FILE_NAME = "hanoi_weather_data.csv"

load_data.py

@@ -7,12 +7,14 @@ from typing import Union
 
 import pandas as pd
 
+from global_config import DATA_FILE_NAME
+
 
 def cut_headers(stream: Union[io.StringIO, io.FileIO, io.TextIOWrapper]):
     headers = []
     exist_startheader = False
-
-    while line := stream.readline():
+    p99_header = 123 # should there be header
+    while (line := stream.readline()) and p99_header > 0:
         if "-BEGIN HEADER-" in line and not exist_startheader:
             exist_startheader = True
             headers.append(str(line.strip()))
@@ -23,8 +25,9 @@ def cut_headers(stream: Union[io.StringIO, io.FileIO, io.TextIOWrapper]):
             # print(f"{line = }")
             headers.append(str(line.strip()))
         elif not exist_startheader:
-            continue
-    if not exist_startheader:
+            p99_header -= 1
+            # continue
+    if not exist_startheader or p99_header == 0:
         # If no headers are found, reset the stream and return
         stream.seek(0)
         return "", stream
@@ -53,27 +56,40 @@ def load_data(
     return data
 
 
-def date_combine(data: pd.DataFrame):
+def date_combine(data: pd.DataFrame, verbose=False):
     data["Date"] = pd.to_datetime(data[["Year", "Month", "Day"]])
     data = data.drop(columns=["Year", "Month", "Day"])
+    if verbose:
+        print("b4 set index", data.columns)
     data = data.set_index("Date")
+    if verbose:
+        print("after set index", data.columns)
+
     return data
 
 
-def validate_view(data: pd.DataFrame):
+def validate_view(data: pd.DataFrame, verbose=False):
+    if verbose:
         print(f"Shape with N/A values: {data.shape}")
         isn = data.isnull().sum().iloc[0]
         print(f"Count of N/A values: {isn}")
     data = data.dropna()  # Drop the rows with missing values
+    if verbose:
         isn = data.isnull().sum().iloc[0]
         print(f"Count after drop: {isn}")
         print(data.describe())
     return data
 
 
-# if __name__ == "__main__":
-def load_this_data(filename):
-    d = load_data(filename)
-    d = date_combine(d)
-    d = validate_view(d)
+def load_this_data(filename=DATA_FILE_NAME, verbose=False, *a):
+    d = load_data(filename, *a)
+    d = date_combine(d, verbose)
+    d = validate_view(d, verbose)
+    if verbose:
+        print(d.head())
+        print(d.columns)
     return d
+
+
+if __name__ == "__main__":
+    load_this_data(verbose=True)

pytorch_port.py

@@ -0,0 +1,183 @@
+import os
+
+from sklearn.model_selection import train_test_split
+import torch
+
+device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
+
+from sklearn.preprocessing import MinMaxScaler
+
+import numpy as np
+from load_data import load_this_data
+
+a = load_this_data()
+
+# time series forecasting
+
+
+scaler = MinMaxScaler()
+print(a.columns)
+scaled_data = scaler.fit_transform(a)
+
+
+def create_sequences(data, seq_length):
+    "see this? return[0] is the big x and return[1] is the big y. and h(x) = y. and h(x) is the model."
+    xs, ys, dates = [], [], []
+    for i in range(len(data) - seq_length):
+        xs.append(data[i : i + seq_length])
+        ys.append(data[i + seq_length])
+        dates.append(i + seq_length)  # index of target
+    return np.array(xs), np.array(ys), np.array(dates)
+
+
+seq_length = 20
+
+X, y, y_indices = create_sequences(scaled_data, seq_length)
+
+X_train, X_test, y_train, y_test = train_test_split(
+    X, y, test_size=0.2, random_state=42
+)
+
+from torch.utils.data import TensorDataset, DataLoader
+
+X_train_tensor = torch.tensor(X_train, dtype=torch.float32)
+y_train_tensor = torch.tensor(y_train, dtype=torch.float32)
+
+train_dataset = TensorDataset(X_train_tensor, y_train_tensor)
+train_loader = DataLoader(train_dataset, batch_size=32, shuffle=True)
+
+
+import torch.nn as nn
+
+
+class LSTMModel(nn.Module):
+    def __init__(self, input_size, hidden_size, num_layers, output_size):
+        super(LSTMModel, self).__init__()
+        self.hidden_size = hidden_size
+        self.num_layers = num_layers
+        self.lstm = nn.LSTM(input_size, hidden_size, num_layers, batch_first=True)
+        self.fc = nn.Linear(hidden_size, output_size)
+
+    def forward(self, x):
+        h0 = torch.zeros(self.num_layers, x.size(0), self.hidden_size).to(device)
+        c0 = torch.zeros(self.num_layers, x.size(0), self.hidden_size).to(device)
+
+        out, _ = self.lstm(x, (h0, c0))
+        out = self.fc(out[:, -1, :])  # Take the output from the last time step
+        return out
+
+
+modfile = "data/model.pth"
+yesmod = os.path.isfile(modfile)
+if yesmod:
+    model = torch.load(modfile)
+else:
+    model = LSTMModel(
+        input_size=X_train.shape[2],
+        hidden_size=64,
+        num_layers=2,
+        output_size=y_train.shape[1],
+    )
+model.to(device)
+
+train_losses = []
+val_losses = []
+
+criterion = nn.MSELoss()
+optimizer = torch.optim.Adam(model.parameters(), lr=0.001)
+
+if not yesmod:
+
+    for epoch in range(100):
+        model.train()
+        running_train_loss = 0.0
+
+        # inputs, targets
+        for batch_X, batch_y in train_loader:
+            batch_X = batch_X.to(device)
+            batch_y = batch_y.to(device)
+
+            outputs = model(batch_X)
+            loss = criterion(outputs, batch_y)
+
+            optimizer.zero_grad()
+            loss.backward()
+            optimizer.step()
+
+            running_train_loss += loss.item() * batch_X.size(0)
+        if (epoch + 1) % 10 == 0:
+            print(f"Epoch [{epoch+1}/100], Loss: {loss.item():.4f}")
+
+        epoch_train_loss = running_train_loss / len(train_loader.dataset)
+        train_losses.append(epoch_train_loss)
+
+    torch.save(model, modfile)
+    # yesmod = True
+
+model.eval()
+
+X_test_tensor = torch.tensor(X_test, dtype=torch.float32).to(device)
+y_test_tensor = torch.tensor(y_test, dtype=torch.float32).to(device)
+val_dataset = TensorDataset(X_test_tensor, y_test_tensor)
+val_loader = DataLoader(val_dataset, batch_size=64, shuffle=False)
+
+# with torch.no_grad():
+#     running_val_loss = 0.0
+#     predictions = model(X_test_tensor)
+#     test_loss = criterion(predictions, y_test_tensor)
+#     print(f"Test Loss: {test_loss.item():.4f}")
+
+running_val_loss = 0.0
+predictions = []
+with torch.no_grad():
+    for inputs, targets in val_loader:
+        inputs, targets = inputs.to(device), targets.to(device)
+        outputs = model(inputs)
+        loss = criterion(outputs, targets)
+        running_val_loss += loss.item() * inputs.size(0)
+        predictions.append(outputs.cpu().numpy())  # Collect all batches
+
+epoch_val_loss = running_val_loss / len(val_loader.dataset)  # len(val_dataset)
+val_losses.append(epoch_val_loss)
+
+predictions = np.concatenate(predictions, axis=0)
+import matplotlib.pyplot as plt
+
+y_test_cpu = y_test_tensor.cpu().numpy()
+
+predictions_rescaled = scaler.inverse_transform(predictions)
+y_test_rescaled = scaler.inverse_transform(y_test_cpu)
+
+
+fig, (ax1, ax2) = plt.subplots(
+    1, 2, figsize=(14, 5), gridspec_kw={"width_ratios": [1, 3]}
+)
+
+# Loss curve on the left
+ax1.plot(train_losses, label="Train Loss")
+ax1.plot(val_losses, label="Val Loss")
+ax1.set_title("Loss Curve")
+ax1.set_xlabel("Epoch")
+ax1.set_ylabel("Loss")
+ax1.legend()
+ax1.grid(True)
+
+# Prediction vs Actual on the right
+
+# 2024 - 25
+_, X_test_idx, _, y_test_idx = train_test_split(
+    X, y_indices, test_size=0.2, random_state=42
+)
+dates_test = a["Date"].iloc[y_test_idx].values  # Actual datetime values
+
+mask = (dates_test >= "2024-01-01") & (dates_test <= "2025-12-31")
+# Plot the results
+ax2.plot(dates_test[mask], y_test_rescaled[mask], label="Actual")
+ax2.plot(dates_test[mask], predictions_rescaled[mask], label="Predicted")
+ax2.set_title("Actual vs. Predicted")
+ax2.set_xlabel("Sample")
+ax2.set_ylabel("Value")
+ax2.legend()
+
+plt.tight_layout()
+plt.show()