MAE_ATMO/torch_MAE_1d_final_40.ipynb

958 lines
182 KiB
Plaintext
Raw Normal View History

2024-11-21 14:02:33 +08:00
{
"cells": [
{
"cell_type": "code",
"execution_count": 25,
"id": "6603a8fc-d9da-4037-b845-d9c38bae4ce4",
"metadata": {},
"outputs": [],
"source": [
"import os\n",
"import torch\n",
"import torch.nn as nn\n",
"import torch.nn.functional as F\n",
"import torch.optim as optim\n",
"from torch.utils.data import DataLoader, Dataset, random_split\n",
"from PIL import Image\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import cv2\n",
"import pandas as pd"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "c28cc123-71be-47ff-b78f-3a4d5592df39",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Maximum pixel value in the dataset: 107.49169921875\n"
]
}
],
"source": [
"\n",
"max_pixel_value = 107.49169921875\n",
"\n",
"print(f\"Maximum pixel value in the dataset: {max_pixel_value}\")"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "dbfe80ce-4394-449c-a9a4-22ed15b2b8f2",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"checkpoint before Generator is OK\n"
]
}
],
"source": [
"class NO2Dataset(Dataset):\n",
" \n",
" def __init__(self, image_dir, mask_dir):\n",
" \n",
" self.image_dir = image_dir\n",
" self.mask_dir = mask_dir\n",
" self.image_filenames = [f for f in os.listdir(image_dir) if f.endswith('.npy')] # 仅加载 .npy 文件\n",
" self.mask_filenames = [f for f in os.listdir(mask_dir) if f.endswith('.jpg')] # 仅加载 .jpg 文件\n",
" \n",
" def __len__(self):\n",
" \n",
" return len(self.image_filenames)\n",
" \n",
" def __getitem__(self, idx):\n",
" \n",
" image_path = os.path.join(self.image_dir, self.image_filenames[idx])\n",
" mask_idx = np.random.choice(self.mask_filenames)\n",
" mask_path = os.path.join(self.mask_dir, mask_idx)\n",
"\n",
" # 加载图像数据 (.npy 文件)\n",
" image = np.load(image_path).astype(np.float32)[:,:,:1] / max_pixel_value # 形状为 (96, 96, 1)\n",
"\n",
" # 加载掩码数据 (.jpg 文件)\n",
" mask = np.array(Image.open(mask_path).convert('L')).astype(np.float32)\n",
"\n",
" # 将掩码数据中非0值设为10值保持不变\n",
" mask = np.where(mask != 0, 1.0, 0.0)\n",
"\n",
" # 保持掩码数据形状为 (96, 96, 1)\n",
" mask = mask[:, :, np.newaxis] # 将形状调整为 (96, 96, 1)\n",
"\n",
" # 应用掩码\n",
" masked_image = image.copy()\n",
" masked_image[:, :, 0] = image[:, :, 0] * mask.squeeze() # 遮盖NO2数据\n",
"\n",
" # cGAN的输入和目标\n",
" X = masked_image[:, :, :1] # 形状为 (96, 96, 8)\n",
" y = image[:, :, 0:1] # 目标输出为NO2数据形状为 (96, 96, 1)\n",
"\n",
" # 转换形状为 (channels, height, width)\n",
" X = np.transpose(X, (2, 0, 1)) # 转换为 (1, 96, 96)\n",
" y = np.transpose(y, (2, 0, 1)) # 转换为 (1, 96, 96)\n",
" mask = np.transpose(mask, (2, 0, 1)) # 转换为 (1, 96, 96)\n",
"\n",
" return torch.tensor(X, dtype=torch.float32), torch.tensor(y, dtype=torch.float32), torch.tensor(mask, dtype=torch.float32)\n",
"\n",
"# 实例化数据集和数据加载器\n",
"image_dir = './out_mat/96/train/'\n",
"mask_dir = './out_mat/96/mask/40/'\n",
"\n",
"print(f\"checkpoint before Generator is OK\")"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "41da7319-9795-441d-bde8-8cf390365099",
"metadata": {},
"outputs": [],
"source": [
"dataset = NO2Dataset(image_dir, mask_dir)\n",
"dataloader = DataLoader(dataset, batch_size=64, shuffle=True, num_workers=8)\n",
"val_set = NO2Dataset('./out_mat/96/valid/', mask_dir)\n",
"val_loader = DataLoader(val_set, batch_size=64, shuffle=False, num_workers=4)\n",
"test_set = NO2Dataset('./out_mat/96/test/', mask_dir)\n",
"test_loader = DataLoader(test_set, batch_size=64, shuffle=False, num_workers=4)"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "70797703-1619-4be7-b965-5506b3d1e775",
"metadata": {},
"outputs": [],
"source": [
"# 可视化特定特征的函数\n",
"def visualize_feature(input_feature,masked_feature, output_feature, title):\n",
" plt.figure(figsize=(12, 6))\n",
" plt.subplot(1, 3, 1)\n",
" plt.imshow(input_feature[0].cpu().numpy(), cmap='RdYlGn_r')\n",
" plt.title(title + \" Input\")\n",
" plt.subplot(1, 3, 2)\n",
" plt.imshow(masked_feature[0].cpu().numpy(), cmap='RdYlGn_r')\n",
" plt.title(title + \" Masked\")\n",
" plt.subplot(1, 3, 3)\n",
" plt.imshow(output_feature[0].detach().cpu().numpy(), cmap='RdYlGn_r')\n",
" plt.title(title + \" Recovery\")\n",
" plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "645114e8-65a4-4867-b3fe-23395288e855",
"metadata": {},
"outputs": [],
"source": [
"class Conv(nn.Sequential):\n",
" def __init__(self, in_channels, out_channels, kernel_size=3, dilation=1, stride=1, bias=False):\n",
" super(Conv, self).__init__(\n",
" nn.Conv2d(in_channels, out_channels, kernel_size=kernel_size, bias=bias,\n",
" dilation=dilation, stride=stride, padding=((stride - 1) + dilation * (kernel_size - 1)) // 2)\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "2af52d0e-b785-4a84-838c-6fcfe2568722",
"metadata": {},
"outputs": [],
"source": [
"class ConvBNReLU(nn.Sequential):\n",
" def __init__(self, in_channels, out_channels, kernel_size=3, dilation=1, stride=1, norm_layer=nn.BatchNorm2d,\n",
" bias=False):\n",
" super(ConvBNReLU, self).__init__(\n",
" nn.Conv2d(in_channels, out_channels, kernel_size=kernel_size, bias=bias,\n",
" dilation=dilation, stride=stride, padding=((stride - 1) + dilation * (kernel_size - 1)) // 2),\n",
" norm_layer(out_channels),\n",
" nn.ReLU()\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "31ecf247-e98b-4977-a145-782914a042bd",
"metadata": {},
"outputs": [],
"source": [
"class SeparableBNReLU(nn.Sequential):\n",
" def __init__(self, in_channels, out_channels, kernel_size=3, stride=1, dilation=1, norm_layer=nn.BatchNorm2d):\n",
" super(SeparableBNReLU, self).__init__(\n",
" nn.Conv2d(in_channels, in_channels, kernel_size=kernel_size, stride=stride, dilation=dilation,\n",
" padding=((stride - 1) + dilation * (kernel_size - 1)) // 2, groups=in_channels, bias=False),\n",
" # 分离卷积,仅调整空间信息\n",
" norm_layer(in_channels), # 对输入通道进行归一化\n",
" nn.Conv2d(in_channels, out_channels, kernel_size=1, bias=False), # 这里进行升维操作\n",
" nn.ReLU6()\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "7827bee2-74f7-4e47-b8c6-e41d5670e8b9",
"metadata": {},
"outputs": [],
"source": [
"class ResidualBlock(nn.Module):\n",
" def __init__(self, in_channels, out_channels, stride=1, downsample=None):\n",
" super(ResidualBlock, self).__init__()\n",
" self.conv1 = nn.Conv2d(in_channels, out_channels, kernel_size=3, stride=stride, padding=1, bias=False)\n",
" self.bn1 = nn.BatchNorm2d(out_channels)\n",
" self.relu = nn.ReLU(inplace=True)\n",
" self.conv2 = nn.Conv2d(out_channels, out_channels, kernel_size=3, padding=1, bias=False)\n",
" self.bn2 = nn.BatchNorm2d(out_channels)\n",
"\n",
" # 如果输入和输出通道不一致,进行降采样操作\n",
" self.downsample = downsample\n",
" if in_channels != out_channels or stride != 1:\n",
" self.downsample = nn.Sequential(\n",
" nn.Conv2d(in_channels, out_channels, kernel_size=1, stride=stride, bias=False),\n",
" nn.BatchNorm2d(out_channels)\n",
" )\n",
"\n",
" def forward(self, x):\n",
" identity = x\n",
" if self.downsample is not None:\n",
" identity = self.downsample(x)\n",
"\n",
" out = self.conv1(x)\n",
" out = self.bn1(out)\n",
" out = self.relu(out)\n",
"\n",
" out = self.conv2(out)\n",
" out = self.bn2(out)\n",
"\n",
" out += identity\n",
" out = self.relu(out)\n",
" return out\n"
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "7853bf62-02f5-4917-b950-6fdfe467df4a",
"metadata": {},
"outputs": [],
"source": [
"class Mlp(nn.Module):\n",
" def __init__(self, in_features, hidden_features=None, out_features=None, act_layer=nn.ReLU6, drop=0.):\n",
" super().__init__()\n",
" out_features = out_features or in_features\n",
" hidden_features = hidden_features or in_features\n",
" self.fc1 = nn.Conv2d(in_features, hidden_features, 1, 1, 0, bias=True)\n",
"\n",
" self.act = act_layer()\n",
" self.fc2 = nn.Conv2d(hidden_features, out_features, 1, 1, 0, bias=True)\n",
" self.drop = nn.Dropout(drop, inplace=True)\n",
"\n",
" def forward(self, x):\n",
" x = self.fc1(x)\n",
" x = self.act(x)\n",
" x = self.drop(x)\n",
" x = self.fc2(x)\n",
" x = self.drop(x)\n",
" return x"
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "e2375881-a11b-47a7-8f56-2eadb25010b0",
"metadata": {},
"outputs": [],
"source": [
"class MultiHeadAttentionBlock(nn.Module):\n",
" def __init__(self, embed_dim, num_heads, dropout=0.1):\n",
" super(MultiHeadAttentionBlock, self).__init__()\n",
" self.attention = nn.MultiheadAttention(embed_dim, num_heads, dropout=dropout)\n",
" self.norm = nn.LayerNorm(embed_dim)\n",
" self.dropout = nn.Dropout(dropout)\n",
"\n",
" def forward(self, x):\n",
" # (B, C, H, W) -> (HW, B, C) for MultiheadAttention compatibility\n",
" B, C, H, W = x.shape\n",
" x = x.view(B, C, H * W).permute(2, 0, 1) # (B, C, H, W) -> (HW, B, C)\n",
"\n",
" # Apply multihead attention\n",
" attn_output, _ = self.attention(x, x, x)\n",
"\n",
" # Apply normalization and dropout\n",
" attn_output = self.norm(attn_output)\n",
" attn_output = self.dropout(attn_output)\n",
"\n",
" # Reshape back to (B, C, H, W)\n",
" attn_output = attn_output.permute(1, 2, 0).view(B, C, H, W)\n",
"\n",
" return attn_output"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "82a15d3d-2f8d-42ec-9146-87c8a4abe384",
"metadata": {},
"outputs": [],
"source": [
"class SpatialAttentionBlock(nn.Module):\n",
" def __init__(self):\n",
" super(SpatialAttentionBlock, self).__init__()\n",
" self.conv = nn.Conv2d(2, 1, kernel_size=7, padding=3, bias=False)\n",
"\n",
" def forward(self, x): #(B, 64, H, W)\n",
" avg_out = torch.mean(x, dim=1, keepdim=True) #(B, 1, H, W)\n",
" max_out, _ = torch.max(x, dim=1, keepdim=True)#(B, 1, H, W)\n",
" out = torch.cat([avg_out, max_out], dim=1)#(B, 2, H, W)\n",
" out = torch.sigmoid(self.conv(out))#(B, 1, H, W)\n",
" return x * out #(B, C, H, W)"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "497bb9f1-1ac5-4d7f-a930-0ea222b9d1d9",
"metadata": {},
"outputs": [],
"source": [
"class DecoderAttentionBlock(nn.Module):\n",
" def __init__(self, in_channels):\n",
" super(DecoderAttentionBlock, self).__init__()\n",
" self.conv1 = nn.Conv2d(in_channels, in_channels // 2, kernel_size=1)\n",
" self.conv2 = nn.Conv2d(in_channels // 2, in_channels, kernel_size=1)\n",
" self.spatial_attention = SpatialAttentionBlock()\n",
"\n",
" def forward(self, x):\n",
" # 通道注意力\n",
" b, c, h, w = x.size()\n",
" avg_pool = F.adaptive_avg_pool2d(x, 1)\n",
" max_pool = F.adaptive_max_pool2d(x, 1)\n",
"\n",
" avg_out = self.conv1(avg_pool)\n",
" max_out = self.conv1(max_pool)\n",
"\n",
" out = avg_out + max_out\n",
" out = torch.sigmoid(self.conv2(out))\n",
"\n",
" # 添加空间注意力\n",
" out = x * out\n",
" out = self.spatial_attention(out)\n",
" return out"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "15b9d453-d8d9-43b8-aca2-904735fb3a99",
"metadata": {},
"outputs": [],
"source": [
"class SEBlock(nn.Module):\n",
" def __init__(self, in_channels, reduced_dim):\n",
" super(SEBlock, self).__init__()\n",
" self.se = nn.Sequential(\n",
" nn.AdaptiveAvgPool2d(1), # 全局平均池化\n",
" nn.Conv2d(in_channels, reduced_dim, kernel_size=1),\n",
" nn.ReLU(),\n",
" nn.Conv2d(reduced_dim, in_channels, kernel_size=1),\n",
" nn.Sigmoid() # 使用Sigmoid是因为我们要对通道进行权重归一化\n",
" )\n",
"\n",
" def forward(self, x):\n",
" return x * self.se(x)"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "6379adb7-8a87-4dd8-a695-4013a7b37830",
"metadata": {
"tags": []
},
"outputs": [],
"source": [
"# 定义Masked Autoencoder模型\n",
"class MaskedAutoencoder(nn.Module):\n",
" def __init__(self):\n",
" super(MaskedAutoencoder, self).__init__()\n",
" self.encoder = nn.Sequential(\n",
" Conv(1, 32, kernel_size=3, stride=2),\n",
" \n",
" nn.ReLU(),\n",
" \n",
" SEBlock(32,32),\n",
" \n",
" ConvBNReLU(32, 64, kernel_size=3, stride=2),\n",
" \n",
" ResidualBlock(64,64),\n",
" \n",
" SeparableBNReLU(64, 128, kernel_size=3, stride=2),\n",
" \n",
" MultiHeadAttentionBlock(embed_dim=128, num_heads=4),\n",
" \n",
" SEBlock(128, 128)\n",
" )\n",
" self.mlp = Mlp(in_features=128, hidden_features=256, out_features=128, act_layer=nn.ReLU6, drop=0.1)\n",
" self.decoder = nn.Sequential(\n",
" nn.ConvTranspose2d(128, 32, kernel_size=3, stride=2, padding=1, output_padding=1),\n",
" nn.ReLU(),\n",
" \n",
" DecoderAttentionBlock(32),\n",
" nn.ConvTranspose2d(32, 16, kernel_size=3, stride=2, padding=1, output_padding=1),\n",
" nn.ReLU(),\n",
" \n",
" DecoderAttentionBlock(16),\n",
" nn.ReLU(),\n",
" \n",
" nn.ConvTranspose2d(16, 1, kernel_size=3, stride=2, padding=1, output_padding=1), # 修改为 output_padding=1\n",
" nn.Sigmoid()\n",
" )\n",
"\n",
" def forward(self, x):\n",
" encoded = self.encoder(x)\n",
" decoded = self.decoder(encoded)\n",
" return decoded\n",
"\n",
"# 实例化模型、损失函数和优化器\n",
"model = MaskedAutoencoder()\n",
"criterion = nn.MSELoss()\n",
"optimizer = torch.optim.Adam(model.parameters(), lr=1e-4)"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "e9c804e0-6f5c-40a7-aba7-a03a496cf427",
"metadata": {},
"outputs": [],
"source": [
"def masked_mse_loss(preds, target, mask):\n",
" loss = (preds - target) ** 2\n",
" loss = loss.mean(dim=-1) # 对每个像素点求平均\n",
" loss = (loss * mask).sum() / mask.sum() # 只计算被mask的像素点的损失\n",
" return loss"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "404a8bfb-4976-4cce-b989-c5e401bce0d7",
"metadata": {},
"outputs": [],
"source": [
"# 训练函数\n",
"def train_epoch(model, device, data_loader, criterion, optimizer):\n",
" model.train()\n",
" running_loss = 0.0\n",
" for batch_idx, (X, y, mask) in enumerate(data_loader):\n",
" X, y, mask = X.to(device), y.to(device), mask.to(device)\n",
" optimizer.zero_grad()\n",
" reconstructed = model(X)\n",
" loss = masked_mse_loss(reconstructed, y, mask)\n",
" # loss = criterion(reconstructed, y)\n",
" loss.backward()\n",
" optimizer.step()\n",
" running_loss += loss.item()\n",
" return running_loss / (batch_idx + 1)"
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "94457c6b-4c6e-4aff-946d-fe4c670bfe16",
"metadata": {},
"outputs": [],
"source": [
"# 评估函数\n",
"def evaluate(model, device, data_loader, criterion):\n",
" model.eval()\n",
" running_loss = 0.0\n",
" with torch.no_grad():\n",
" for batch_idx, (X, y, mask) in enumerate(data_loader):\n",
" X, y, mask = X.to(device), y.to(device), mask.to(device)\n",
" reconstructed = model(X)\n",
" if batch_idx == 8:\n",
" rand_ind = np.random.randint(0, len(y))\n",
" # visualize_feature(y[rand_ind], X[rand_ind], reconstructed[rand_ind], title='NO_2')\n",
" loss = masked_mse_loss(reconstructed, y, mask)\n",
" running_loss += loss.item()\n",
" return running_loss / (batch_idx + 1)"
]
},
{
"cell_type": "code",
"execution_count": 19,
"id": "296ba6bd-2239-4948-b278-7edcb29bfd14",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"cuda\n"
]
}
],
"source": [
"# 数据准备\n",
"device = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n",
"print(device)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "743d1000-561e-4444-8b49-88346c14f28b",
"metadata": {},
"outputs": [],
"source": [
"model = model.to(device)\n",
"\n",
"num_epochs = 100\n",
"train_losses = list()\n",
"val_losses = list()\n",
"for epoch in range(num_epochs):\n",
" train_loss = train_epoch(model, device, dataloader, criterion, optimizer)\n",
" train_losses.append(train_loss)\n",
" val_loss = evaluate(model, device, val_loader, criterion)\n",
" val_losses.append(val_loss)\n",
" print(f'Epoch {epoch+1}, Train Loss: {train_loss}, Val Loss: {val_loss}')\n",
"\n",
"# 测试模型\n",
"test_loss = evaluate(model, device, test_loader, criterion)\n",
"print(f'Test Loss: {test_loss}')"
]
},
{
"cell_type": "code",
"execution_count": 38,
"id": "cdc0d608-6f0a-43dc-8cc1-8acf68215d18",
"metadata": {},
"outputs": [
{
"ename": "NameError",
"evalue": "name 'train_losses' is not defined",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mNameError\u001b[0m Traceback (most recent call last)",
"Cell \u001b[0;32mIn[38], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m tr_ind \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mlist\u001b[39m(\u001b[38;5;28mrange\u001b[39m(\u001b[38;5;28mlen\u001b[39m(\u001b[43mtrain_losses\u001b[49m)))\n\u001b[1;32m 2\u001b[0m val_ind \u001b[38;5;241m=\u001b[39m \u001b[38;5;28mlist\u001b[39m(\u001b[38;5;28mrange\u001b[39m(\u001b[38;5;28mlen\u001b[39m(val_losses)))\n\u001b[1;32m 3\u001b[0m plt\u001b[38;5;241m.\u001b[39mplot(train_losses[\u001b[38;5;241m1\u001b[39m:], label\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mtrain_loss\u001b[39m\u001b[38;5;124m'\u001b[39m)\n",
"\u001b[0;31mNameError\u001b[0m: name 'train_losses' is not defined"
]
}
],
"source": [
"tr_ind = list(range(len(train_losses)))\n",
"val_ind = list(range(len(val_losses)))\n",
"plt.plot(train_losses[1:], label='train_loss')\n",
"plt.plot(val_losses[1:], label='val_loss')\n",
"plt.legend(loc='best')"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "849b1706-1a98-4571-989f-da06d949c843",
"metadata": {},
"outputs": [],
"source": [
"torch.save(model, './models/MAE/final_40.pt')"
]
},
{
"cell_type": "code",
"execution_count": 20,
"id": "40a803b2-4891-4d47-ab61-cf88db8007a0",
"metadata": {},
"outputs": [],
"source": [
"model = torch.load('./models/MAE/final_40.pt')"
]
},
{
"cell_type": "code",
"execution_count": 21,
"id": "a8467686-0655-4056-8e01-56299eb89d7c",
"metadata": {},
"outputs": [],
"source": [
"from sklearn.metrics import mean_squared_error, r2_score, mean_absolute_percentage_error, mean_absolute_error"
]
},
{
"cell_type": "code",
"execution_count": 22,
"id": "016c3045-0312-462f-82ae-7272944ed92d",
"metadata": {},
"outputs": [],
"source": [
"def cal_ioa(y_true, y_pred):\n",
" # 计算平均值\n",
" mean_observed = np.mean(y_true)\n",
" mean_predicted = np.mean(y_pred)\n",
"\n",
" # 计算IoA\n",
" numerator = np.sum((y_true - y_pred) ** 2)\n",
" denominator = np.sum((np.abs(y_true - mean_observed) + np.abs(y_pred - mean_predicted)) ** 2)\n",
" IoA = 1 - (numerator / denominator)\n",
"\n",
" return IoA"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "dae7427e-548e-4276-a4ea-bc9b279d44e8",
"metadata": {},
"outputs": [],
"source": [
"eva_list = list()\n",
"device = 'cpu'\n",
"model = model.to(device)\n",
"with torch.no_grad():\n",
" for batch_idx, (X, y, mask) in enumerate(test_loader):\n",
" X, y, mask = X.to(device), y.to(device), mask.to(device)\n",
" mask_rev = (torch.squeeze(mask, dim=1)==0) * 1 # mask取反获得修复区域\n",
" reconstructed = model(X)\n",
" rev_data = y * max_pixel_value\n",
" rev_recon = reconstructed * max_pixel_value\n",
" # todo: 这里需要只评估修补出来的模块\n",
" data_label = torch.squeeze(rev_data, dim=1) * mask_rev\n",
" data_label = data_label[mask_rev==1]\n",
" recon_no2 = torch.squeeze(rev_recon, dim=1) * mask_rev\n",
" recon_no2 = recon_no2[mask_rev==1]\n",
" mae = mean_absolute_error(data_label, recon_no2)\n",
" rmse = np.sqrt(mean_squared_error(data_label, recon_no2))\n",
" mape = mean_absolute_percentage_error(data_label, recon_no2)\n",
" r2 = r2_score(data_label, recon_no2)\n",
" ioa = cal_ioa(data_label.detach().numpy(), recon_no2.detach().numpy())\n",
" r = np.corrcoef(data_label, recon_no2)[0, 1]\n",
" eva_list.append([mae, rmse, mape, r2, ioa, r])"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "d5b01834-ca18-4ec3-bc9d-64382d0fab34",
"metadata": {},
"outputs": [],
"source": [
"pd.DataFrame(eva_list, columns=['mae', 'rmse', 'mape', 'r2', 'ioa', 'r']).describe()"
]
},
{
"cell_type": "code",
"execution_count": 23,
"id": "0887481a-764e-4fd5-9580-45aa813a4391",
"metadata": {},
"outputs": [],
"source": [
"eva_list_frame = list()\n",
"device = 'cpu'\n",
"model = model.to(device)\n",
"best_mape = 1\n",
"best_img = None\n",
"best_mask = None\n",
"best_recov = None\n",
"with torch.no_grad():\n",
" for batch_idx, (X, y, mask) in enumerate(test_loader):\n",
" X, y, mask = X.to(device), y.to(device), mask.to(device)\n",
" mask_rev = (torch.squeeze(mask, dim=1)==0) * 1 # mask取反获得修复区域\n",
" reconstructed = model(X)\n",
" rev_data = y * max_pixel_value\n",
" rev_recon = reconstructed * max_pixel_value\n",
" # todo: 这里需要只评估修补出来的模块\n",
" for i, sample in enumerate(rev_data):\n",
" used_mask = mask_rev[i]\n",
" data_label = sample[0] * used_mask\n",
" recon_no2 = rev_recon[i][0] * used_mask\n",
" data_label = data_label[used_mask==1]\n",
" recon_no2 = recon_no2[used_mask==1]\n",
" mae = mean_absolute_error(data_label, recon_no2)\n",
" rmse = np.sqrt(mean_squared_error(data_label, recon_no2))\n",
" mape = mean_absolute_percentage_error(data_label, recon_no2)\n",
" r2 = r2_score(data_label, recon_no2)\n",
" ioa = cal_ioa(data_label.detach().numpy(), recon_no2.detach().numpy())\n",
" r = np.corrcoef(data_label, recon_no2)[0, 1]\n",
" eva_list_frame.append([mae, rmse, mape, r2, ioa, r])\n",
" if mape < best_mape:\n",
" best_recov = rev_recon[i][0].numpy()\n",
" best_mask = used_mask.numpy()\n",
" best_img = sample[0].numpy()\n",
" best_mape = mape"
]
},
{
"cell_type": "code",
"execution_count": 26,
"id": "f7355895-ffde-458f-b4e6-b8afd95ea663",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th>mae</th>\n",
" <th>rmse</th>\n",
" <th>mape</th>\n",
" <th>r2</th>\n",
" <th>ioa</th>\n",
" <th>r</th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th>count</th>\n",
" <td>4739.000000</td>\n",
" <td>4739.000000</td>\n",
" <td>4739.000000</td>\n",
" <td>4739.000000</td>\n",
" <td>4739.000000</td>\n",
" <td>4739.000000</td>\n",
" </tr>\n",
" <tr>\n",
" <th>mean</th>\n",
" <td>1.540401</td>\n",
" <td>2.199879</td>\n",
" <td>0.195554</td>\n",
" <td>0.585799</td>\n",
" <td>0.848016</td>\n",
" <td>0.778401</td>\n",
" </tr>\n",
" <tr>\n",
" <th>std</th>\n",
" <td>0.647315</td>\n",
" <td>0.909418</td>\n",
" <td>0.092239</td>\n",
" <td>0.213993</td>\n",
" <td>0.106987</td>\n",
" <td>0.127430</td>\n",
" </tr>\n",
" <tr>\n",
" <th>min</th>\n",
" <td>0.462070</td>\n",
" <td>0.593854</td>\n",
" <td>0.068942</td>\n",
" <td>-0.551587</td>\n",
" <td>0.218504</td>\n",
" <td>0.145717</td>\n",
" </tr>\n",
" <tr>\n",
" <th>25%</th>\n",
" <td>1.021385</td>\n",
" <td>1.472757</td>\n",
" <td>0.144170</td>\n",
" <td>0.460184</td>\n",
" <td>0.805011</td>\n",
" <td>0.711952</td>\n",
" </tr>\n",
" <tr>\n",
" <th>50%</th>\n",
" <td>1.367836</td>\n",
" <td>2.056208</td>\n",
" <td>0.176119</td>\n",
" <td>0.624006</td>\n",
" <td>0.876770</td>\n",
" <td>0.805993</td>\n",
" </tr>\n",
" <tr>\n",
" <th>75%</th>\n",
" <td>1.975825</td>\n",
" <td>2.792777</td>\n",
" <td>0.217419</td>\n",
" <td>0.745375</td>\n",
" <td>0.923944</td>\n",
" <td>0.871612</td>\n",
" </tr>\n",
" <tr>\n",
" <th>max</th>\n",
" <td>5.186517</td>\n",
" <td>9.158884</td>\n",
" <td>0.960081</td>\n",
" <td>0.968376</td>\n",
" <td>0.992196</td>\n",
" <td>0.985054</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"</div>"
],
"text/plain": [
" mae rmse mape r2 ioa \\\n",
"count 4739.000000 4739.000000 4739.000000 4739.000000 4739.000000 \n",
"mean 1.540401 2.199879 0.195554 0.585799 0.848016 \n",
"std 0.647315 0.909418 0.092239 0.213993 0.106987 \n",
"min 0.462070 0.593854 0.068942 -0.551587 0.218504 \n",
"25% 1.021385 1.472757 0.144170 0.460184 0.805011 \n",
"50% 1.367836 2.056208 0.176119 0.624006 0.876770 \n",
"75% 1.975825 2.792777 0.217419 0.745375 0.923944 \n",
"max 5.186517 9.158884 0.960081 0.968376 0.992196 \n",
"\n",
" r \n",
"count 4739.000000 \n",
"mean 0.778401 \n",
"std 0.127430 \n",
"min 0.145717 \n",
"25% 0.711952 \n",
"50% 0.805993 \n",
"75% 0.871612 \n",
"max 0.985054 "
]
},
"execution_count": 26,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"pd.DataFrame(eva_list_frame, columns=['mae', 'rmse', 'mape', 'r2', 'ioa', 'r']).describe()"
]
},
{
"cell_type": "code",
"execution_count": 27,
"id": "ee712c69-2b57-4ac6-9c7d-d73ba0d1ecca",
"metadata": {},
"outputs": [],
"source": [
"# 可视化特定特征的函数\n",
"def visualize_rst(input_feature,masked_feature, recov_region, output_feature, title):\n",
" plt.figure(figsize=(12, 6))\n",
" plt.subplot(1, 4, 1)\n",
" plt.imshow(input_feature, cmap='RdYlGn_r')\n",
" plt.gca().axis('off') # 获取当前坐标轴并关闭\n",
" plt.subplot(1, 4, 2)\n",
" plt.imshow(masked_feature, cmap='gray')\n",
" plt.gca().axis('off') # 获取当前坐标轴并关闭\n",
" plt.subplot(1, 4, 3)\n",
" plt.imshow(recov_region, cmap='RdYlGn_r')\n",
" plt.gca().axis('off') # 获取当前坐标轴并关闭\n",
" plt.subplot(1, 4, 4)\n",
" plt.imshow(output_feature, cmap='RdYlGn_r')\n",
" plt.gca().axis('off') # 获取当前坐标轴并关闭\n",
" plt.savefig('./figures/result/40_samples.png', bbox_inches='tight')"
]
},
{
"cell_type": "code",
"execution_count": 28,
"id": "e89feac0-d03a-4686-8a38-722e6a54a96f",
"metadata": {},
"outputs": [],
"source": [
"best_mask_cp = np.where(best_mask == 0, np.nan, best_mask)"
]
},
{
"cell_type": "code",
"execution_count": 34,
"id": "a87835de-836b-411b-b4cc-68e98b6638f4",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[30.36338043, 30.67309189, 30.94369125, ..., 11.77855492,\n",
" 11.96412849, 11.9506712 ],\n",
" [30.04488182, 30.25416946, 30.87792015, ..., 11.70056629,\n",
" 12.05164337, 11.96099949],\n",
" [29.82366371, 30.49637985, 30.7125721 , ..., 11.49174881,\n",
" 11.77280235, 11.96125317],\n",
" ...,\n",
" [ 8.4842186 , 9.02253723, 8.97320557, ..., 5.35319471,\n",
" 5.15942717, 5.25348282],\n",
" [ 8.59376144, 8.57794476, 8.91248322, ..., 5.41437721,\n",
" 5.41615629, 5.49798965],\n",
" [ 8.4524231 , 8.80022049, 8.73760223, ..., 5.64806128,\n",
" 5.53445244, 5.61840296]])"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"best_recov * (1-best_mask) + best_recov*best_mask"
]
},
{
"cell_type": "code",
"execution_count": 37,
"id": "95e3882b-9962-4aab-be80-4240f326ef51",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA7YAAADeCAYAAAAJtZwyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy88F64QAAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOy9eZwlVXk+/pxzqqrv9PSsDAOjwDCssoiKSIyKESRuwd1I3HFJRIPGGPfoV1QkrlERFTQmGBf4uSVG8yWuURRRvoiCuIAwAwyyzzBLT/ftWs75/fGe9yxVdbtvbzPdw33m03PvrVvLqbr11nmX531fYYwxGGCAAQYYYIABBhhggAEGGGCARQq5pwcwwAADDDDAAAMMMMAAAwwwwACzwcCwHWCAAQYYYIABBhhggAEGGGBRY2DYDjDAAAMMMMAAAwwwwAADDLCoMTBsBxhggAEGGGCAAQYYYIABBljUGBi2AwwwwAADDDDAAAMMMMAAAyxqDAzbAQYYYIABBhhggAEGGGCAARY1BobtAAMMMMAAAwwwwAADDDDAAIsaA8N2gAEGGGCAAQYYYIABBhhggEWNgWE7wAADDDDAAAMMMMAAAwwwwKJG0u+K5edfCAAwhQYqAygBABDSvnYSIJX0WUogS4HOkN24BLQBuhOA1jClpn3tKmAKDTNRwlQGKDSM1jA7C+idOaqJErrQ0IVGOVHBaANdVBBSIB1OIaRAMpxASAHZSSBSCUgBoSSNJZUQqYLoKIihBHJZBjGkgOElQKKA4Q6NVQogSYCVyyCSDBCS/lTi3wOA0YAu6X1V+mUAwNt1ltM+J8boO1HzHRjt96O1Xy4lIJN4nyrzy1VCrzIBOiMoTBnvFhqFzjFW7ICGRqlzGqYpoY1GobsodYWxYgKl0ahDGwNdWy6DsYffZSqBhEAiFaQQbvlIOoIlyQiEkFAiifahjYYxGrnuQhsaX2VK3DexE92yQK5LlFqj1BW0MShNhdxe48oYaEO3Hb2aaJwTpYG9faCNgRQCSgCFNig0UFT0yvsK0bH3MX9vb2ekyu6jon1PBPsIx1RoYDSvUGnj9v2Vv/hc4/ouJNx22219rffpT38a55xzDkztmj30oQ/Fpz/9aaxbtw4AoJTCPvvsgzRNsXXrVoyPj7t1P/ShD+FjH/sYXvrSl+Kd73wnlFIAgE6ng9WrV0PKmfvWvva1r+GMM87A6OgoAODYY4/F17/+dRx++OF972NsbAxbt251n4eGhrDPPvvMalw/+tGPcNZZZ+Huu+/Gli1bUFXVjPc1wJ5F/d4foB2Xrz8SQDylDS+habUzLJEOp0hHUiRDCtlIhnQkhRpKIIYTDJ97aes+Jz72HAglaE6387tIFc3tiSIdw87f8s/Pb92H/ubf0KCkBDqZ102koAHyeymBREEkCYzVV+Thb53z6zTAAG3Q338tvSntXMGCVFb0XmsYbXrr12UFMzoOU1QwOwuYokK1tQvTLWHGSuiiIl261FBDCmooQTlWoBwvIaQgOYPX5xl8TDWkIFMFuSwDUgk5nJJePZzS+6FQHqXXr8NzCuQMSeL17kTRHwBICSGtzm3IHoDWfh+MLIVgnVhI0r9D3V3YcQCkqxsNFF0+KXod6Nf2HPY+/bpvw1akdOMZe3QhBRm3StCNmEqILPUbsMGojb3BNN28WkKghNEm3rckwYWGncgUlDbQhYZMJVIlyPhdYn9Qa0SrVNH6oVHLV68yQEcEY69ovKW9abs5rZulNM7RMZhkAmJ4iTcyZSAouiQh0SWMvSmcEBZd2iYf8+uyEUsrBgKXxEas0fQ5yYILHhjWocDa7yQk3cyg/bNgZKqDQucokTuhAwApEqRSIpElENrTdp+lqVBq2g/fvHVoe2P1+r4yJXLdhRIJjNDQ0Kh0CSkkEpnBGA0hJCQAISQEJCQEpBD+VQh3HH6vhABAAlAZg7GCvlciFEYSNHtRnVBUBuhWBhOVvy3oHOizcrcKbcPgfZEwx9+7Y2qDvNLIq+aDbCHjL//yL/ta74477mhV7G+88UaceeaZyDK6X/fZZx+cc845OOKII/ChD30IP/zhD926mzdvBgBceuml+N3vfueWn3jiiXj3u9+NFStWzOJMZo8f/OAHeP/734/SPhOOO+44nHPOOdh3331nvM+HPvSh+Ld/+zdce+21+Md//EfceeedczXcAQZYkBhZ6vRv8PTKOq6QAqayDmoAQhZ2mUGimvMIQ987DjGcuPkblX3VmvQPNlgnc0JJAUjl12G9xDm0SakWaebmZqGCeXuAAXYHdpBzVo8XLnAkpPQBpGQK/ToBxPAQRFnBpAqmqKDHaF/yAUNIpCQjt9CkB2uDRJtYD5cCybA9hl1eFWRQylR5ZQkgg7MQQLeEBiC0JvonG64A6ddaA3lBn9nw7WS0f2fMCkBbx5JKyKAtcmfkGgCQgR4i7bWJrkXidevQ4GXdOnxlPXqgX++1+nXfhq07yxr4BnOeHtl7onIGbguEkjCo/LGk9dQG+3PR4cDD5K5e3ahtQxW4Jcoy9iABzjtmspQMZaOtoW231zqIuLILQvpJUJPR6wSKt0GwD/YkRVFgO0GrxK/jDODAEG71BNnot30VQkIKGXmAjNEkqAL25vYPBRYgaYR9L6FN5fYde5W8wGk07wUDDW00BDQAEvrKlDCQUCa+1aSQMEbWlsW/ne4RLWEBCAWh/h0v19o4waH5QtjvDADhvEB1oeLbqNDemxUeI/QgAfRZTXbvLSD87Gc/a12+ZMkSZ6wCwMTEROt6o6OjuPrqq93ntWvX4vbbb8eaNWtw3XXXte7/zjvvjAy8TqeDLVu2QCmF4eHhWUVIGVpr7NixA9u3b+97m5tvvhlXXnkl8jx3+7j33nvddRBCYOnSpS7S3A9WrFiBE044AUoprFmzBrt27cKuXbug9fQf0AMMsBjwkN9eDwD49bFHumV1m9NUBkgpClQVmubwordMLHnPf896XPIvLow+6/99XWNwIqkpxXZONps/CugSYv0bZj2OAQaYDKYkK8F0K6AyEKmEUQYCsUE5pX6dkJEJ2IBRqomt2EmgKUQGdAVMt4z0a6MNyaOiZWzwKgC6CoJZgNejZWD1VGQkC2gyUlm/1sZHWxPvmHIPBm1ipS2EXU9ICSN1rJgxQkYk69Zs1Aay3OqoCiO2A/06WB6/Lkb9uv+IbScjT6n2VOTQo0QjCqxArf3ZSOndueGVUQICkoSo8mcupAA6iR9guI30wgclIIYSek3tWHLthNJoA1H5X8lYOrMEKMJcVn5sSeJfs9R6iTwdIfQQmUDgkHRiirKpCZGUaPUkSUmCxVFg/gx4AzcU0EAAKlM6b5Ex5Lnh10qXDaFzQ4FEphIksukZkkIgMcZSFdgLJd13UqiGIEhBdAmmWYhayraBFTwjUSCP9gkAduSRp8p9Z4zz6gQ/Ifj2C8HbjRW0jQoEuDKWLqH5VjWB4JLwhejlOapTLEJkStIMsIiRJAn+5m/+Bqeeeqpb9vWvfx0XXXTRlHTM7du34z3veQ9WrlyJX/3qV30d77rrrsOrXvUqHHXUUXjrW9+K/fbbbzbDB0DR4Te84Q0YGRnpe5tbb70VRVG4z3/4wx/w2te+Fp0O0ahWr16Nt7zlLTjqqKOmPZ5DDjkEH/vYx3DTTTfh3HPPxc033zztfQwwwGJC1pFIEj/nsJKstYHUBqaiP11qFGPEyFq6OwfINEimP2apV2yTDEjtfM70RV1Ovc8BBpgleurXTLsHrP6MQJ+eXL+WIxnMkPKMB2By/Vp6vVra701RTaraUAohRYBNoSldoKNpX6xfJ5YxwfTjMAWAx2/Pz8DTPUL92mlpMjBZQt1Y1fRrpibrQIbdq30+caR2oF/bZXuPft1/xDZRgLa0XqYFhFFTRuSNCaKaUrorZ9q8NCrwBikbka3I4xTTlq2XlSnIHWUNW0s3rgwEEG1jglC2gKUks1dLaqAMrhp7mJIEKAs3foPKR3jdztjrK2Oh4QcN4I3UulHLE2nl6Q/Rd4B/reXphl4kFjoniPZzLyiRQAptvUX+d/MCZBoepzZaBH0
"text/plain": [
"<Figure size 1200x600 with 4 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"visualize_rst(best_img, best_mask, best_recov*best_mask_cp, best_img * (1-best_mask) + best_recov*best_mask, '')"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c32a4b1f-9e2d-46cd-b117-9857dc840c7c",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.16"
}
},
"nbformat": 4,
"nbformat_minor": 5
}