diff --git a/Baseline_models/ssb_models/models/agcrn/__init__.py b/Baseline_models/ssb_models/models/agcrn/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/Baseline_models/ssb_models/models/agcrn/__init__.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/agcrn/__init__.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/Baseline_models/ssb_models/models/agcrn/layer.py b/Baseline_models/ssb_models/models/agcrn/layer.py new file mode 100644 index 0000000000000000000000000000000000000000..0bc31d0cf79e86ec61d3445601f41fc448f98113 --- /dev/null +++ b/Baseline_models/ssb_models/models/agcrn/layer.py @@ -0,0 +1,62 @@ +import torch +import torch.nn as nn +import torch.nn.functional as F + + +class AVWGCN(nn.Module): + def __init__(self, dim_in, dim_out, cheb_k, embed_dim): + super(AVWGCN, self).__init__() + self.cheb_k = cheb_k + self.weights_pool = nn.Parameter( + torch.zeros(embed_dim, cheb_k, dim_in, dim_out) + ) + self.bias_pool = nn.Parameter(torch.zeros(embed_dim, dim_out)) + nn.init.xavier_uniform_(self.weights_pool) + nn.init.xavier_uniform_(self.bias_pool) + + def forward(self, x, node_embeddings): + # x shaped[B, N, C], node_embeddings shaped [N, D] -> supports shaped [N, N] + # output shape [B, N, C] + node_num = node_embeddings.shape[0] + supports = F.softmax( + F.relu(torch.mm(node_embeddings, node_embeddings.transpose(0, 1))), dim=1 + ) + support_set = [torch.eye(node_num).to(supports.device), supports] + # default cheb_k = 3 + for k in range(2, self.cheb_k): + support_set.append( + torch.matmul(2 * supports, support_set[-1]) - support_set[-2] + ) + supports = torch.stack(support_set, dim=0) + weights = torch.einsum( + "nd,dkio->nkio", node_embeddings, self.weights_pool + ) # N, cheb_k, dim_in, dim_out + bias = torch.matmul(node_embeddings, self.bias_pool) # N, dim_out + x_g = torch.einsum("knm,bmc->bknc", supports, x) # B, cheb_k, N, dim_in + x_g = x_g.permute(0, 2, 1, 3) # B, N, cheb_k, dim_in + x_gconv = torch.einsum("bnki,nkio->bno", x_g, weights) + bias # b, N, dim_out + return x_gconv + + +class AGCRNCell(nn.Module): + def __init__(self, node_num, dim_in, dim_out, cheb_k, embed_dim): + super(AGCRNCell, self).__init__() + self.node_num = node_num + self.hidden_dim = dim_out + self.gate = AVWGCN(dim_in + self.hidden_dim, 2 * dim_out, cheb_k, embed_dim) + self.update = AVWGCN(dim_in + self.hidden_dim, dim_out, cheb_k, embed_dim) + + def forward(self, x, state, node_embeddings): + # x: B, num_nodes, input_dim + # state: B, num_nodes, hidden_dim + state = state.to(x.device) + input_and_state = torch.cat((x, state), dim=-1) + z_r = torch.sigmoid(self.gate(input_and_state, node_embeddings)) + z, r = torch.split(z_r, self.hidden_dim, dim=-1) + candidate = torch.cat((x, z * state), dim=-1) + hc = torch.tanh(self.update(candidate, node_embeddings)) + h = r * state + (1 - r) * hc + return h + + def init_hidden_state(self, batch_size): + return torch.zeros(batch_size, self.node_num, self.hidden_dim) diff --git "a/Baseline_models/ssb_models/models/agcrn/layer.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/agcrn/layer.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/Baseline_models/ssb_models/models/agcrn/net.py b/Baseline_models/ssb_models/models/agcrn/net.py new file mode 100644 index 0000000000000000000000000000000000000000..e27b3cd2a1ec2b97d48a08e7fbc8e6bfcfb4479d --- /dev/null +++ b/Baseline_models/ssb_models/models/agcrn/net.py @@ -0,0 +1,113 @@ +import torch +import torch.nn as nn + +from Baseline_models.ssb_models.models.agcrn.layer import AGCRNCell + + +class AVWDCRNN(nn.Module): + def __init__( + self, + node_num: int, + dim_in: int, + dim_out: int, + cheb_k: int, + embed_dim: int, + num_layers: int, + dropout: int, + ): + super(AVWDCRNN, self).__init__() + assert num_layers >= 1, "At least one DCRNN layer in the Encoder." + self.node_num = node_num + self.input_dim = dim_in + self.num_layers = num_layers + self.dropout = nn.Dropout(dropout) + self.dcrnn_cells = nn.ModuleList() + self.dcrnn_cells.append(AGCRNCell(node_num, dim_in, dim_out, cheb_k, embed_dim)) + for _ in range(1, num_layers): + self.dcrnn_cells.append( + AGCRNCell(node_num, dim_out, dim_out, cheb_k, embed_dim) + ) + + def forward(self, x, init_state, node_embeddings): + # shape of x: (B, T, N, D) + # shape of init_state: (num_layers, B, N, hidden_dim) + assert x.shape[2] == self.node_num and x.shape[3] == self.input_dim + seq_length = x.shape[1] + current_inputs = x + output_hidden = [] + for i in range(self.num_layers): + state = init_state[i] + inner_states = [] + for t in range(seq_length): + state = self.dropout( + self.dcrnn_cells[i]( + current_inputs[:, t, :, :], state, node_embeddings + ) + ) + inner_states.append(state) + output_hidden.append(state) + current_inputs = torch.stack(inner_states, dim=1) + # current_inputs: the outputs of last layer: (B, T, N, hidden_dim) + # output_hidden: the last state for each layer: (num_layers, B, N, hidden_dim) + # last_state: (B, N, hidden_dim) + return current_inputs, output_hidden + + def init_hidden(self, batch_size): + init_states = [] + for i in range(self.num_layers): + init_states.append(self.dcrnn_cells[i].init_hidden_state(batch_size)) + return torch.stack(init_states, dim=0) # (num_layers, B, N, hidden_dim) + + +class AGCRN(nn.Module): + def __init__( + self, + num_nodes: int, + input_dim: int, + rnn_units: int, + output_dim: int, + window_size: int, + dropout: int, + num_layers: int, + embedding_dim: int, + cheb_k: int, + ): + super(AGCRN, self).__init__() + self.num_nodes = num_nodes + self.input_dim = input_dim + self.hidden_dim = rnn_units + self.output_dim = output_dim + self.window_size = window_size + self.num_layers = num_layers + + self.node_embeddings = nn.Parameter( + torch.randn(self.num_nodes, embedding_dim), requires_grad=True + ) + + self.encoder = AVWDCRNN( + num_nodes, input_dim, rnn_units, cheb_k, embedding_dim, num_layers, dropout + ) + + # predictor + self.end_conv = nn.Conv2d( + 1, + window_size * self.output_dim, + kernel_size=(1, self.hidden_dim), + bias=True, + ) + + def forward(self, x): + # source: B, T_1, N, D + # target: B, T_2, N, D + # supports = F.softmax(F.relu(torch.mm(self.nodevec1, self.nodevec1.transpose(0,1))), dim=1) + + init_state = self.encoder.init_hidden(x.shape[0]) + output, _ = self.encoder(x, init_state, self.node_embeddings) # B, T, N, hidden + output = output[:, -1:, :, :] # B, 1, N, hidden + + # CNN based predictor + output = self.end_conv(output) # B, T*C, N, 1 + output = output.squeeze(-1).reshape( + -1, self.window_size, self.output_dim, self.num_nodes + ) + return output.permute(0, 1, 3, 2) # B, T, N, C diff --git "a/Baseline_models/ssb_models/models/agcrn/net.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/agcrn/net.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/Baseline_models/ssb_models/models/mtgnn/__init__.py b/Baseline_models/ssb_models/models/mtgnn/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/Baseline_models/ssb_models/models/mtgnn/__init__.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/mtgnn/__init__.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/Baseline_models/ssb_models/models/mtgnn/layer.py b/Baseline_models/ssb_models/models/mtgnn/layer.py new file mode 100644 index 0000000000000000000000000000000000000000..b0bfe54bc98ae99a7611e7f31c1b9d77f932b9b4 --- /dev/null +++ b/Baseline_models/ssb_models/models/mtgnn/layer.py @@ -0,0 +1,258 @@ +from __future__ import division + +from typing import Optional + +import torch +import torch.nn as nn +from torch.nn import init +import numbers +import torch.nn.functional as F + + +class NConv(nn.Module): + def __init__(self): + super(NConv, self).__init__() + + def forward(self, x, A): + x = torch.einsum("ncwl,vw->ncvl", (x, A)) + return x.contiguous() + + +class DyNConv(nn.Module): + def __init__(self): + super(DyNConv, self).__init__() + + def forward(self, x, A): + x = torch.einsum("ncvl,nvwl->ncwl", (x, A)) + return x.contiguous() + + +class Linear(nn.Module): + def __init__(self, c_in, c_out, bias=True): + super(Linear, self).__init__() + self.mlp = torch.nn.Conv2d( + c_in, c_out, kernel_size=(1, 1), padding=(0, 0), stride=(1, 1), bias=bias + ) + + def forward(self, x): + return self.mlp(x) + + +class Propagation(nn.Module): + def __init__(self, c_in, c_out, gdep, dropout, alpha): + super(Propagation, self).__init__() + self.nconv = NConv() + self.mlp = Linear(c_in, c_out) + self.gdep = gdep + self.dropout = dropout + self.alpha = alpha + + def forward(self, x, adj): + adj = adj + torch.eye(adj.size(0)).to(x.device) + d = adj.sum(1) + h = x + dv = d + a = adj / dv.view(-1, 1) + for i in range(self.gdep): + h = self.alpha * x + (1 - self.alpha) * self.nconv(h, a) + ho = self.mlp(h) + return ho + + +class MixPropagation(nn.Module): + def __init__(self, c_in, c_out, gdep, dropout, alpha): + super(MixPropagation, self).__init__() + self.nconv = NConv() + self.mlp = Linear((gdep + 1) * c_in, c_out) + self.gdep = gdep + self.dropout = dropout + self.alpha = alpha + + def forward(self, x, adj): + adj = adj + torch.eye(adj.size(0)).to(x.device) + d = adj.sum(1) + h = x + out = [h] + a = adj / d.view(-1, 1) + for i in range(self.gdep): + h = self.alpha * x + (1 - self.alpha) * self.nconv(h, a) + out.append(h) + ho = torch.cat(out, dim=1) + ho = self.mlp(ho) + return ho + + +class DyMixPropagation(nn.Module): + def __init__(self, c_in, c_out, gdep, dropout, alpha): + super(DyMixPropagation, self).__init__() + self.nconv = DyNConv() + self.mlp1 = Linear((gdep + 1) * c_in, c_out) + self.mlp2 = Linear((gdep + 1) * c_in, c_out) + + self.gdep = gdep + self.dropout = dropout + self.alpha = alpha + self.lin1 = Linear(c_in, c_in) + self.lin2 = Linear(c_in, c_in) + + def forward(self, x): + x1 = torch.tanh(self.lin1(x)) + x2 = torch.tanh(self.lin2(x)) + adj = self.nconv(x1.transpose(2, 1), x2) + adj0 = torch.softmax(adj, dim=2) + adj1 = torch.softmax(adj.transpose(2, 1), dim=2) + + h = x + out = [h] + for i in range(self.gdep): + h = self.alpha * x + (1 - self.alpha) * self.nconv(h, adj0) + out.append(h) + ho = torch.cat(out, dim=1) + ho1 = self.mlp1(ho) + + h = x + out = [h] + for i in range(self.gdep): + h = self.alpha * x + (1 - self.alpha) * self.nconv(h, adj1) + out.append(h) + ho = torch.cat(out, dim=1) + ho2 = self.mlp2(ho) + + return ho1 + ho2 + + +class Dilated1DConv(nn.Module): + def __init__(self, cin, cout, dilation_factor=2): + super(Dilated1DConv, self).__init__() + self.kernel_set = [2, 3, 6, 7] + self.tconv = nn.Conv2d(cin, cout, (1, 7), dilation=(1, dilation_factor)) + + def forward(self, input): + x = self.tconv(input) + return x + + +class DilatedInception(nn.Module): + def __init__(self, cin, cout, dilation_factor=2): + super(DilatedInception, self).__init__() + self.tconv = nn.ModuleList() + self.kernel_set = [2, 3, 6, 7] + cout = int(cout / len(self.kernel_set)) + for kern in self.kernel_set: + self.tconv.append( + nn.Conv2d(cin, cout, (1, kern), dilation=(1, dilation_factor)) + ) + + def forward(self, input): + x = [] + for i in range(len(self.kernel_set)): + x.append(self.tconv[i](input)) + for i in range(len(self.kernel_set)): + x[i] = x[i][..., -x[-1].size(3) :] + x = torch.cat(x, dim=1) + return x + + +class GraphConstructor(nn.Module): + def __init__( + self, + nnodes: int, + k: int, + dim: int, + alpha=3, + static_feat: Optional[torch.Tensor] = None, + ): + super(GraphConstructor, self).__init__() + self.nnodes = nnodes + if static_feat is not None: + xd = static_feat.shape[1] + self.lin1 = nn.Linear(xd, dim) + self.lin2 = nn.Linear(xd, dim) + else: + self.emb1 = nn.Embedding(nnodes, dim) + self.emb2 = nn.Embedding(nnodes, dim) + self.lin1 = nn.Linear(dim, dim) + self.lin2 = nn.Linear(dim, dim) + + self.k = k + self.dim = dim + self.alpha = alpha + self.static_feat = static_feat + + def forward(self, idx): + if self.static_feat is None: + nodevec1 = self.emb1(idx) + nodevec2 = self.emb2(idx) + else: + nodevec1 = self.static_feat[idx, :] + nodevec2 = nodevec1 + + nodevec1 = torch.tanh(self.alpha * self.lin1(nodevec1)) + nodevec2 = torch.tanh(self.alpha * self.lin2(nodevec2)) + + a = torch.mm(nodevec1, nodevec2.transpose(1, 0)) - torch.mm( + nodevec2, nodevec1.transpose(1, 0) + ) + adj = F.relu(torch.tanh(self.alpha * a)) + mask = torch.zeros( + idx.size(0), idx.size(0), device=idx.device, dtype=torch.float + ) + s1, t1 = (adj + torch.rand_like(adj) * 0.01).topk(self.k, 1) + mask.scatter_(1, t1, s1.fill_(1)) + adj = adj * mask + return adj + + def fullA(self, idx): + if self.static_feat is None: + nodevec1 = self.emb1(idx) + nodevec2 = self.emb2(idx) + else: + nodevec1 = self.static_feat[idx, :] + nodevec2 = nodevec1 + + nodevec1 = torch.tanh(self.alpha * self.lin1(nodevec1)) + nodevec2 = torch.tanh(self.alpha * self.lin2(nodevec2)) + + a = torch.mm(nodevec1, nodevec2.transpose(1, 0)) - torch.mm( + nodevec2, nodevec1.transpose(1, 0) + ) + adj = F.relu(torch.tanh(self.alpha * a)) + return adj + + +class LayerNorm(nn.Module): + __constants__ = ["normalized_shape", "weight", "bias", "eps", "elementwise_affine"] + + def __init__(self, normalized_shape, eps=1e-5, elementwise_affine=True): + super(LayerNorm, self).__init__() + if isinstance(normalized_shape, numbers.Integral): + normalized_shape = (normalized_shape,) + self.normalized_shape = tuple(normalized_shape) + self.eps = eps + self.elementwise_affine = elementwise_affine + if self.elementwise_affine: + self.weight = nn.Parameter(torch.ones(*normalized_shape)) + self.bias = nn.Parameter(torch.zeros(*normalized_shape)) + else: + self.register_parameter("weight", None) + self.register_parameter("bias", None) + + def forward(self, input, idx): + if self.elementwise_affine: + return F.layer_norm( + input, + tuple(input.shape[1:]), + self.weight[:, idx, :], + self.bias[:, idx, :], + self.eps, + ) + else: + return F.layer_norm( + input, tuple(input.shape[1:]), self.weight, self.bias, self.eps + ) + + def extra_repr(self): + return ( + "{normalized_shape}, eps={eps}, " + "elementwise_affine={elementwise_affine}".format(**self.__dict__) + ) diff --git "a/Baseline_models/ssb_models/models/mtgnn/layer.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/mtgnn/layer.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/Baseline_models/ssb_models/models/mtgnn/net.py b/Baseline_models/ssb_models/models/mtgnn/net.py new file mode 100644 index 0000000000000000000000000000000000000000..58256be18ac4dfec044ec08e3e6b23bc445b6db1 --- /dev/null +++ b/Baseline_models/ssb_models/models/mtgnn/net.py @@ -0,0 +1,259 @@ +from typing import Optional + +import torch +from torch import nn +import torch.nn.functional as F + +import pytorch_lightning as pl + +from Baseline_models.ssb_models.models.mtgnn.layer import ( + GraphConstructor, + DilatedInception, + MixPropagation, + LayerNorm, +) + + +#class MtGNN(pl.LightningModule): +class MtGNN(nn.Module): + def __init__( + self, + gcn_depth: int, + num_nodes: int, + window_size: int, + static_feat=None, + dropout=0.3, + subgraph_size=20, + node_dim=40, + dilation_exponential=1, + conv_channels=16, + residual_channels=16, + skip_channels=16, + end_channels=16, + in_dim=2, + layers=3, + propalpha=0.05, + tanhalpha=3, + layer_norm_affine=True, + out_dim=8 + ): + super(MtGNN, self).__init__() + self.num_nodes = num_nodes + self.dropout = dropout + self.filter_convs = nn.ModuleList() + self.gate_convs = nn.ModuleList() + self.residual_convs = nn.ModuleList() + self.skip_convs = nn.ModuleList() + self.gconv1 = nn.ModuleList() + self.gconv2 = nn.ModuleList() + self.norm = nn.ModuleList() + self.start_conv = nn.Conv2d( + in_channels=in_dim, out_channels=residual_channels, kernel_size=(1, 1) + ) + self.gc = GraphConstructor( + num_nodes, + subgraph_size, + node_dim, + alpha=tanhalpha, + static_feat=static_feat, + ) + + self.seq_length = window_size + kernel_size = 7 + if dilation_exponential > 1: + self.receptive_field = int( + 1 + + (kernel_size - 1) + * (dilation_exponential**layers - 1) + / (dilation_exponential - 1) + ) + else: + self.receptive_field = layers * (kernel_size - 1) + 1 + + for i in range(1): + if dilation_exponential > 1: + rf_size_i = int( + 1 + + i + * (kernel_size - 1) + * (dilation_exponential**layers - 1) + / (dilation_exponential - 1) + ) + else: + rf_size_i = i * layers * (kernel_size - 1) + 1 + new_dilation = 1 + for j in range(1, layers + 1): + if dilation_exponential > 1: + rf_size_j = int( + rf_size_i + + (kernel_size - 1) + * (dilation_exponential**j - 1) + / (dilation_exponential - 1) + ) + else: + rf_size_j = rf_size_i + j * (kernel_size - 1) + + self.filter_convs.append( + DilatedInception( + residual_channels, conv_channels, dilation_factor=new_dilation + ) + ) + self.gate_convs.append( + DilatedInception( + residual_channels, conv_channels, dilation_factor=new_dilation + ) + ) + self.residual_convs.append( + nn.Conv2d( + in_channels=conv_channels, + out_channels=residual_channels, + kernel_size=(1, 1), + ) + ) + if self.seq_length > self.receptive_field: + self.skip_convs.append( + nn.Conv2d( + in_channels=conv_channels, + out_channels=skip_channels, + kernel_size=(1, self.seq_length - rf_size_j + 1), + ) + ) + else: + self.skip_convs.append( + nn.Conv2d( + in_channels=conv_channels, + out_channels=skip_channels, + kernel_size=(1, self.receptive_field - rf_size_j + 1), + ) + ) + + self.gconv1.append( + MixPropagation( + conv_channels, + residual_channels, + gcn_depth, + dropout, + propalpha, + ) + ) + self.gconv2.append( + MixPropagation( + conv_channels, + residual_channels, + gcn_depth, + dropout, + propalpha, + ) + ) + + if self.seq_length > self.receptive_field: + self.norm.append( + LayerNorm( + ( + residual_channels, + num_nodes, + self.seq_length - rf_size_j + 1, + ), + elementwise_affine=layer_norm_affine, + ) + ) + else: + self.norm.append( + LayerNorm( + ( + residual_channels, + num_nodes, + self.receptive_field - rf_size_j + 1, + ), + elementwise_affine=layer_norm_affine, + ) + ) + + new_dilation *= dilation_exponential + + self.layers = layers + self.end_conv_1 = nn.Conv2d( + in_channels=skip_channels, + out_channels=end_channels, + kernel_size=(1, 1), + bias=True, + ) + self.end_conv_2 = nn.Conv2d( + in_channels=end_channels, + out_channels=out_dim, + kernel_size=(1, 1), + bias=True, + ) + if self.seq_length > self.receptive_field: + self.skip0 = nn.Conv2d( + in_channels=in_dim, + out_channels=skip_channels, + kernel_size=(1, self.seq_length), + bias=True, + ) + self.skipE = nn.Conv2d( + in_channels=residual_channels, + out_channels=skip_channels, + kernel_size=(1, self.seq_length - self.receptive_field + 1), + bias=True, + ) + + else: + self.skip0 = nn.Conv2d( + in_channels=in_dim, + out_channels=skip_channels, + kernel_size=(1, self.receptive_field), + bias=True, + ) + self.skipE = nn.Conv2d( + in_channels=residual_channels, + out_channels=skip_channels, + kernel_size=(1, 1), + bias=True, + ) + + self.idx = nn.Parameter(torch.arange(self.num_nodes), requires_grad=False) + + def forward(self, input: torch.Tensor, idx=None): + input = input.transpose(1, 3) + seq_len = input.size(3) + assert ( + seq_len == self.seq_length + ), "input sequence length not equal to preset sequence length" + + if self.seq_length < self.receptive_field: + input = nn.functional.pad( + input, (self.receptive_field - self.seq_length, 0, 0, 0) + ) + + if idx is None: + adp = self.gc(self.idx) + else: + adp = self.gc(idx) + + x = self.start_conv(input) + skip = self.skip0(F.dropout(input, self.dropout, training=self.training)) + for i in range(self.layers): + residual = x + filter = self.filter_convs[i](x) + filter = torch.tanh(filter) + gate = self.gate_convs[i](x) + gate = torch.sigmoid(gate) + x = filter * gate + x = F.dropout(x, self.dropout, training=self.training) + s = x + s = self.skip_convs[i](s) + skip = s + skip + x = self.gconv1[i](x, adp) + self.gconv2[i](x, adp.transpose(1, 0)) + + x = x + residual[:, :, :, -x.size(3) :] + if idx is None: + x = self.norm[i](x, self.idx) + else: + x = self.norm[i](x, idx) + + skip = self.skipE(x) + skip + x = F.relu(skip) + x = F.relu(self.end_conv_1(x)) + x = self.end_conv_2(x) + return x diff --git "a/Baseline_models/ssb_models/models/mtgnn/net.py\357\200\272Zone.Identifier" "b/Baseline_models/ssb_models/models/mtgnn/net.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/ModeConvFast.py b/ModeConvModel/ModeConvFast.py new file mode 100644 index 0000000000000000000000000000000000000000..61888275508f9078cf505a1512e8e03b8e683842 --- /dev/null +++ b/ModeConvModel/ModeConvFast.py @@ -0,0 +1,77 @@ +from typing import Optional + +import numpy as np +import scipy.ndimage +from scipy import linalg as LA +import torch +from torch import Tensor +import torch.nn.functional as F +from torch.nn import Parameter + +from torch_geometric.nn.conv import MessagePassing +from torch_geometric.nn.dense.linear import Linear +from torch_geometric.nn.inits import zeros +from torch_geometric.typing import OptTensor +from torch_geometric.utils import ( + add_self_loops, + get_laplacian, + remove_self_loops, +) + + +class ModeConvFast(MessagePassing): + def __init__(self, in_channels: int, out_channels: int, K: int, num_sensors, + bias: bool = True, **kwargs): + kwargs.setdefault('aggr', 'add') + super().__init__(**kwargs) + + self.in_channels = in_channels + self.out_channels = out_channels + self.K = K + self.num_sensors = num_sensors + self.lins_r = torch.nn.ModuleList([ + Linear(in_channels, out_channels, bias=False, + weight_initializer='glorot') for _ in range(K) + ]) + self.lins_i = torch.nn.ModuleList([ + Linear(in_channels, out_channels, bias=False, + weight_initializer='glorot') for _ in range(K) + ]) + + if bias: + self.bias = Parameter(torch.Tensor(out_channels)) + else: + self.register_parameter('bias', None) + + self.reset_parameters() + + def reset_parameters(self): + for lin in self.lins_i: + lin.reset_parameters() + for lin in self.lins_r: + lin.reset_parameters() + zeros(self.bias) + + def forward(self, x: Tensor, edge_index: Tensor, + weight: OptTensor = None, sim=None): + + weight = weight * sim[..., None] + Tx_1 = self.propagate(edge_index, x=x, weight=weight, size=None) + Tx_1_c = torch.view_as_real(Tx_1) + out_r = self.lins_r[0](Tx_1_c[..., 0]).view((-1, self.out_channels)) + out_i = self.lins_i[0](Tx_1_c[..., 1]).view((-1, self.out_channels)) + out = out_r - out_i + + if self.bias is not None: + out += self.bias + + return out + + def message(self, x_j, weight): + r = weight.reshape(-1, 1, 2)[..., 0] * x_j + i = weight.reshape(-1, 1, 2)[..., 1] * x_j + return torch.complex(r, i) + + def __repr__(self) -> str: + return (f'{self.__class__.__name__}({self.in_channels}, ' + f'{self.out_channels}, K={len(self.lins)}') diff --git "a/ModeConvModel/ModeConvFast.py\357\200\272Zone.Identifier" "b/ModeConvModel/ModeConvFast.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/ModeConvLaplace.py b/ModeConvModel/ModeConvLaplace.py new file mode 100644 index 0000000000000000000000000000000000000000..79d14a43c58e6b631bd1da89c777dbf3523178b3 --- /dev/null +++ b/ModeConvModel/ModeConvLaplace.py @@ -0,0 +1,75 @@ +from typing import Optional + +import numpy as np +from scipy import linalg as LA +import torch +from torch import Tensor +from torch.nn import Parameter + +from torch_geometric.nn.conv import MessagePassing +from torch_geometric.nn.dense.linear import Linear +from torch_geometric.nn.inits import zeros +from torch_geometric.typing import OptTensor +from torch_geometric.utils import ( + add_self_loops, + get_laplacian, + remove_self_loops, +) + + +class ModeConvLaplace(MessagePassing): + def __init__(self, in_channels: int, out_channels: int, K: int, num_sensors, + bias: bool = True, **kwargs): + kwargs.setdefault('aggr', 'add') + super().__init__(**kwargs) + + self.in_channels = in_channels + self.out_channels = out_channels + self.K = K + self.num_sensors = num_sensors + self.lins_r = torch.nn.ModuleList([ + Linear(in_channels, out_channels, bias=False, + weight_initializer='glorot') for _ in range(K) + ]) + self.lins_i = torch.nn.ModuleList([ + Linear(in_channels, out_channels, bias=False, + weight_initializer='glorot') for _ in range(K) + ]) + + if bias: + self.bias = Parameter(torch.Tensor(out_channels)) + else: + self.register_parameter('bias', None) + + self.reset_parameters() + + def reset_parameters(self): + for lin in self.lins_i: + lin.reset_parameters() + for lin in self.lins_r: + lin.reset_parameters() + zeros(self.bias) + + def forward(self, x: Tensor, edge_index: Tensor, + weight: OptTensor = None, sim=None): + + weight = weight * sim[..., None] + Tx_1 = self.propagate(edge_index, x=x, weight=weight, size=None) + Tx_1_c = torch.view_as_real(Tx_1) + out_r = self.lins_r[0](Tx_1_c[..., 0]).view((-1, self.out_channels)) + out_i = self.lins_i[0](Tx_1_c[..., 1]).view((-1, self.out_channels)) + out = out_r - out_i + + if self.bias is not None: + out += self.bias + + return out + + def message(self, x_j, weight): + r = weight.reshape(-1, 1, 2)[..., 0] * x_j + i = weight.reshape(-1, 1, 2)[..., 1] * x_j + return torch.complex(r, i) + + def __repr__(self) -> str: + return (f'{self.__class__.__name__}({self.in_channels}, ' + f'{self.out_channels}, K={len(self.lins)}') diff --git "a/ModeConvModel/ModeConvLaplace.py\357\200\272Zone.Identifier" "b/ModeConvModel/ModeConvLaplace.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/dataset.py b/ModeConvModel/dataset.py new file mode 100644 index 0000000000000000000000000000000000000000..7d0095ab3bd465fa61bd90037ab3cb9901c40925 --- /dev/null +++ b/ModeConvModel/dataset.py @@ -0,0 +1,234 @@ +import json + +import networkx +import numpy as np +import pandas as pd +import torch +import tqdm +from torch_geometric.data import Data +from torch_geometric.data import InMemoryDataset +from torch_geometric.utils.convert import from_networkx + + +#create sample +def create_sample( + strain, + temp, + inclination, + cov, + labels, + edge_index, +): + + sample = Data(x=torch.concatenate((inclination, temp), dim=-1), edge_index=edge_index, y=labels) + sample.strain = strain + sample.cov = cov + return sample + + +def create_sample_luxem( + acceleration, + temp, + voltage, + disp, + cov, + labels, + edge_index +): + sample = Data(x=torch.concatenate((voltage, temp), dim=-1), edge_index=edge_index, y=labels) + sample.acc = acceleration + sample.disp = disp + sample.cov = cov + return sample + + +def batch_cov(points): + """https://stackoverflow.com/a/71357620""" + B, N, D = points.size() + mean = points.mean(dim=1).unsqueeze(1) + diffs = (points - mean).reshape(B * N, D) + prods = torch.bmm(diffs.unsqueeze(2), diffs.unsqueeze(1)).reshape(B, N, D, D) + bcov = prods.sum(dim=1) / (N - 1) # Unbiased estimate + return bcov # (B, D, D) + + +class SimulatedSmartBridgeDataset(InMemoryDataset): + def __init__(self, root, mode="train", data="short", transform=None, pre_transform=None, pre_filter=None): + self.mode = mode + self.data = data + super().__init__(root, transform, pre_transform, pre_filter) + self.data, self.slices = torch.load(self.processed_paths[0]) + + @property + def processed_file_names(self): + if self.mode == "train": + return ["train_data.pt"] + elif self.mode == "val": + return ["val_data.pt"] + elif self.mode == "test": + return ["test_data.pt"] + + def process(self): + data = pd.read_csv(f"data/SimulatedSmartBridge/simulated_smart_bridge_{self.mode}.csv") + inclination = data.values[:, :3] + temp = data.values[:, 3:4] + strain = data.values[:, 4:16] + labels = data.values[:, -1] + + num_sensors = 12 + G = networkx.complete_graph(num_sensors) + pyg_graph: Data = from_networkx(G) + + lookback = 5 + tau = 5 + num = lookback * tau + + # write metadata + if self.mode == "train": + stats = {} + stats["inclination"] = {"mean": inclination.mean(), "std": inclination.std(), "min": inclination.min(), "max": inclination.max()} + stats["temp"] = {"mean": temp.mean(), "std": temp.std(), "min": temp.min(), "max": temp.max()} + stats["strain"] = {"mean": strain.mean(), "std": strain.std(), "min": strain.min(), "max": strain.max()} + with open(f'{self.processed_dir}/metadata.json', 'w') as fp: + json.dump(stats, fp, indent=4) + + strain = torch.Tensor(strain.reshape((-1, lookback, tau, num_sensors))) + cov = torch.zeros( + (strain.shape[0] * strain.shape[1], strain.shape[-1], strain.shape[-1])) + strain2 = torch.nn.functional.pad(strain.view((-1, strain.shape[-2], strain.shape[-1])), + (0, 0, 0, 0, 0, 1)) + for i in range(num_sensors): + strain3 = strain2.clone() + strain3[0:-1, :, np.delete(np.arange(num_sensors), i)] = strain3[1:, :, + np.delete(np.arange(num_sensors), i)] + cov[:, i, :] = batch_cov(strain3)[:-1, i] + for i in range(num_sensors): + s = strain2[:-1, :, i].unsqueeze(-1) + s_lag = strain2[1:, :, i].unsqueeze(-1) + s = torch.concatenate((s, s_lag), -1) + test = batch_cov(s) + cov[:, i, i] = test[:, 0, 1] + + cov = cov.view((strain.shape[0], strain.shape[1], strain.shape[-1], strain.shape[-1])) + + cov = (cov - (-370000)) / (370000 - (-370000)) + + strain = strain.reshape((-1, num, num_sensors)) + temp = torch.Tensor(temp.reshape((-1, 1, num, 1))).repeat(1, num_sensors, 1, 1).squeeze(-1) + inclination = torch.Tensor(np.repeat(inclination.reshape((-1, num, 3)), [4,4,4], axis=2).transpose(0,2,1)) + labels = labels[num - 1::num] + data_list = [] + for i in tqdm.tqdm(range(len(labels))): + sample = create_sample( + strain[i], + temp[i], + inclination[i], + cov[i], + labels[i], + pyg_graph.edge_index, + ) + + # Read data into huge `Data` list. + data_list.append(sample) + + + if self.pre_filter is not None: + data_list = [data for data in data_list if self.pre_filter(data)] + + if self.pre_transform is not None: + data_list = [self.pre_transform(data) for data in data_list] + + data, slices = self.collate(data_list) + torch.save((data, slices), self.processed_paths[0]) + + +class LuxemburgDataset(InMemoryDataset): + def __init__(self, root, mode="train", transform=None, pre_transform=None, pre_filter=None): + self.mode = mode + super().__init__(root, transform, pre_transform, pre_filter) + self.data, self.slices = torch.load(self.processed_paths[0]) + + @property + def processed_file_names(self): + if self.mode == "train": + return ["train_data.pt"] + elif self.mode == "val": + return ["val_data.pt"] + elif self.mode == "test": + return ["test_data.pt"] + + def process(self): + data = pd.read_csv(f"data/luxemburg/lux_{self.mode}.csv") + voltage = data.values[:, :4] + temp = data.values[:, 4:5] + disp = data.values[:, 5:13] + acceleration = data.values[:, 13:39] + labels = data.values[:, -1] + + num_sensors = 26 + G = networkx.complete_graph(num_sensors) + pyg_graph: Data = from_networkx(G) + + lookback = 5 + tau = 5 + num = lookback * tau + + # write metadata + if self.mode == "train": + stats = {} + stats["voltage"] = {"mean": list(voltage.mean(0)), "std": list(voltage.std(0)), "min": list(voltage.min(0)), "max": list(voltage.max(0))} + stats["temp"] = {"mean": list(temp.mean(0)), "std": list(temp.std(0)), "min": list(temp.min(0)), "max": list(temp.max(0))} + stats["disp"] = {"mean": list(disp.mean(0)), "std": list(disp.std(0)), "min": list(disp.min(0)), "max": list(disp.max(0))} + stats["acceleration"] = {"mean": acceleration.mean(), "std": acceleration.std(), "min": acceleration.min(), "max": acceleration.max()} + with open(f'{self.processed_dir}/metadata.json', 'w') as fp: + json.dump(stats, fp, indent=4) + + acceleration = torch.Tensor(acceleration.reshape((-1, lookback, tau, num_sensors))) + cov = torch.zeros( + (acceleration.shape[0] * acceleration.shape[1], acceleration.shape[-1], acceleration.shape[-1])) + acceleration2 = torch.nn.functional.pad(acceleration.view((-1, acceleration.shape[-2], acceleration.shape[-1])), + (0, 0, 0, 0, 0, 1)) + for i in range(num_sensors): + acceleration3 = acceleration2.clone() + acceleration3[0:-1, :, np.delete(np.arange(num_sensors), i)] = acceleration3[1:, :, + np.delete(np.arange(num_sensors), i)] + cov[:, i, :] = batch_cov(acceleration3)[:-1, i] + for i in range(num_sensors): + acc = acceleration2[:-1, :, i].unsqueeze(-1) + acc_lag = acceleration2[1:, :, i].unsqueeze(-1) + acc = torch.concatenate((acc, acc_lag), -1) + test = batch_cov(acc) + cov[:, i, i] = test[:, 0, 1] + + cov = cov.view((acceleration.shape[0], acceleration.shape[1], acceleration.shape[-1], acceleration.shape[-1])) + + cov = (cov - (-0.001)) / (0.001 - (-0.001)) + + acceleration = acceleration.reshape((-1, num, num_sensors)) + temp = torch.Tensor(temp.reshape((-1, num, 1))) + voltage = torch.Tensor(voltage.reshape((-1, num, 4))) + disp = torch.Tensor(disp.reshape((-1, num, 8))) + labels = labels[num-1::num] + data_list = [] + for i in tqdm.tqdm(range(len(acceleration))): + sample = create_sample_luxem( + acceleration[i], + temp[i], + voltage[i], + disp[i], + cov[i], + labels[i], + pyg_graph.edge_index, + ) + + # Read data into huge `Data` list. + data_list.append(sample) + + if self.pre_filter is not None: + data_list = [data for data in data_list if self.pre_filter(data)] + + if self.pre_transform is not None: + data_list = [self.pre_transform(data) for data in data_list] + + data, slices = self.collate(data_list) + torch.save((data, slices), self.processed_paths[0]) diff --git "a/ModeConvModel/dataset.py\357\200\272Zone.Identifier" "b/ModeConvModel/dataset.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/LuxAGCRNModelPL.py b/ModeConvModel/models/Luxemburg/LuxAGCRNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..a19636322348cd78e3ddefc39f1bdc9c379e991c --- /dev/null +++ b/ModeConvModel/models/Luxemburg/LuxAGCRNModelPL.py @@ -0,0 +1,185 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics +from Baseline_models.ssb_models.models.agcrn.net import AGCRN + + +class AGCRNEdgeDecoder(torch.nn.Module): + def __init__(self, args): + super().__init__() + self.lin1 = Linear(args.bottleneck, args.hidden_dim // 2) + self.lin2 = Linear(args.hidden_dim // 2, args.hidden_dim) + self.lin3 = Linear(args.hidden_dim, 14 * 25) + + def forward(self, z): + z = self.lin1(z).relu() + z = self.lin2(z).relu() + z = self.lin3(z) + return z + + +class AGCRNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = AGCRN(num_sensors, input_dim=14, rnn_units=args.hidden_dim, output_dim=args.bottleneck, window_size=1, + dropout=0, num_layers=args.num_layer, embedding_dim=2, cheb_k=2) + if args.decoder == "linear": + self.decoder = AGCRNEdgeDecoder(args) + elif args.decoder == "custom": + self.decoder = AGCRN(num_sensors, input_dim=args.bottleneck, rnn_units=args.hidden_dim, output_dim=14*25, window_size=1, + dropout=0, num_layers=args.num_layer, embedding_dim=2, cheb_k=2) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, x): + z = self.encoder(x) + z = self.decoder(z) + return z + + +class LuxAGCRNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = AGCRNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_lux/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.concatenate((stats["voltage"]["min"], stats["temp"]["min"], stats["disp"]["min"], [stats["acceleration"]["min"]])), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.concatenate((stats["voltage"]["max"], stats["temp"]["max"], stats["disp"]["max"], [stats["acceleration"]["max"]])), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.concatenate((stats["voltage"]["mean"], stats["temp"]["mean"], stats["disp"]["mean"], [stats["acceleration"]["mean"]])), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.concatenate((stats["voltage"]["std"], stats["temp"]["std"], stats["disp"]["std"], [stats["acceleration"]["std"]])), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, acceleration, disp = batch.x, batch.edge_index, batch.acc, batch.disp + if batch_idx == 0: + self.epoch_start_time = time.time() + acceleration = acceleration.view((batch.num_graphs, 5, 5, acceleration.shape[-1])).permute(0, 3, 1, 2).reshape((batch.num_graphs, 25, acceleration.shape[-1], 1)) + + num_sensors = acceleration.shape[-2] + disp = disp.view((batch.num_graphs, -1, 1, disp.shape[-1])).repeat((1, 1, num_sensors, 1)) + + x = x.view((batch.num_graphs, -1, 1, x.shape[-1])) + x = x.repeat(1, 1, num_sensors, 1) + x = torch.concatenate((x, disp, acceleration), dim=-1) + + # x = (x - self.mean) / self.std + x = (x - self.min) / (self.max - self.min) + + return x, edge_index + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index).view(train_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + x = x.permute(0, 2, 1, 3) + + losses = torch.nn.MSELoss()(y_hat, x) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y_hat = y_hat.view(val_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + + batchsize = len(val_batch.y) + x = x.permute(0, 2, 1, 3) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.reshape(-1).cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index).view(test_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + y = test_batch.y + + batchsize = len(test_batch.y) + x = x.permute(0, 2, 1, 3) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/Luxemburg/LuxAGCRNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/LuxAGCRNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/LuxChebGNNModelPL.py b/ModeConvModel/models/Luxemburg/LuxChebGNNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..27e3edf9a59c5290b4affeff8587aa2d424e7b32 --- /dev/null +++ b/ModeConvModel/models/Luxemburg/LuxChebGNNModelPL.py @@ -0,0 +1,243 @@ +import json +import os +import pickle +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from sklearn.metrics import confusion_matrix, balanced_accuracy_score, roc_auc_score +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn import ChebConv +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics + + +class ChebGNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ChebConv(14 * 25, args.bottleneck, K=5)]) + self.lin_layers = nn.ModuleList([Linear(14 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ChebConv(14 * 25, args.hidden_dim, K=5)]) + self.lin_layers = nn.ModuleList([Linear(14 * 25, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append(ChebConv(args.hidden_dim, args.hidden_dim, K=5)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ChebConv(args.hidden_dim, args.bottleneck, K=5)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index) + lin(x) + x = x.relu() + return x + + +class ChebEdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ChebConv(args.bottleneck, 14 * 25, K=5, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.conv_layers = nn.ModuleList( + [ChebConv(args.bottleneck, args.hidden_dim, K=5, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append( + ChebConv(args.hidden_dim, args.hidden_dim, K=5, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ChebConv(args.hidden_dim, 14 * 25, K=5, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index) + self.lin_layers[-1](z) + return z.view(-1) + + +class ChebGNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = ChebGNNEncoder(num_sensors, args) + self.decoder = ChebEdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index): + z = self.encoder(x, edge_index) + z = self.decoder(z, edge_index) + return z + + +class LuxChebGNNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = ChebGNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_lux/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.concatenate((np.tile(stats["voltage"]["min"], 25), np.tile(stats["temp"]["min"], 25), np.tile(stats["disp"]["min"], 25), np.tile(stats["acceleration"]["min"], 25))), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.concatenate((np.tile(stats["voltage"]["max"], 25), np.tile(stats["temp"]["max"], 25), np.tile(stats["disp"]["max"], 25), np.tile(stats["acceleration"]["max"], 25))), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.concatenate((np.tile(stats["voltage"]["mean"], 25), np.tile(stats["temp"]["mean"], 25), np.tile(stats["disp"]["mean"], 25), np.tile(stats["acceleration"]["mean"], 25))), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.concatenate((np.tile(stats["voltage"]["std"], 25), np.tile(stats["temp"]["std"], 25), np.tile(stats["disp"]["std"], 25), np.tile(stats["acceleration"]["std"], 25))), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x, edge_index) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, acceleration, disp = batch.x, batch.edge_index, batch.acc, batch.disp + if batch_idx == 0: + self.epoch_start_time = time.time() + acceleration = acceleration.view((batch.num_graphs, 5, 5, acceleration.shape[-1])) + + num_sensors = acceleration.shape[-1] + + e = np.vstack((np.repeat(np.arange(num_sensors), num_sensors), np.tile(np.arange(num_sensors), num_sensors))) + e = e[:, np.where(e[0] != e[1])[0]] + add = np.arange(0, batch.num_graphs) + add = np.repeat(add, e.shape[1]) + add *= num_sensors + e1 = np.tile(e, batch.num_graphs) + e2 = e1 + add + edge_index = torch.tensor(e2, device=self.dev, dtype=torch.int64) + + disp = disp.view((batch.num_graphs, 1, -1, disp.shape[-1])).repeat((1, num_sensors, 1, 1)) + + # reshape and add acceleration to x + x = x.view((batch.num_graphs, 1, -1, x.shape[-1])) + x = x.repeat(1, num_sensors, 1, 1) + x = torch.concatenate((x.view((x.shape[0], x.shape[1], -1)), disp.view((x.shape[0], x.shape[1], -1)), + acceleration.permute(0, 3, 1, 2).reshape(acceleration.shape[0], acceleration.shape[3], + -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + # x = (x - self.mean) / self.std + x = (x - self.min) / (self.max - self.min) + + return x, edge_index + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + + losses = torch.nn.MSELoss()(y_hat, x.view(-1)) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/Luxemburg/LuxChebGNNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/LuxChebGNNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/LuxModeConvFastModelPL.py b/ModeConvModel/models/Luxemburg/LuxModeConvFastModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..1067321368616277f0faa4e493e6acc57a865f14 --- /dev/null +++ b/ModeConvModel/models/Luxemburg/LuxModeConvFastModelPL.py @@ -0,0 +1,325 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy import linalg as LA +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.ModeConvFast import ModeConvFast +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics + + +class GNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvFast(14 * 25, args.bottleneck, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(14 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ModeConvFast(14 * 25, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(14*25, args.hidden_dim)]) + for i in range(args.num_layer-2): + self.conv_layers.append(ModeConvFast(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvFast(args.hidden_dim, args.bottleneck, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index, weight, sim): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index, weight=weight, sim=sim) + lin(x) + x = x.relu() + return x + + +class EdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvFast(args.bottleneck, 14 * 25, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.conv_layers = nn.ModuleList( + [ModeConvFast(args.bottleneck, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append( + ModeConvFast(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvFast(args.hidden_dim, 14 * 25, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index, weight, sim): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index, weight=weight, sim=sim) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index, weight=weight, sim=sim) + self.lin_layers[-1](z) + return z.view(-1) + + +class GNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = GNNEncoder(num_sensors, args) + self.decoder = EdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index, weight, sim): + z = self.encoder(x, edge_index, weight, sim) + z = self.decoder(z, edge_index, weight, sim) + return z + + +class LuxModeConvFastModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = GNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_lux/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.concatenate((np.tile(stats["voltage"]["min"], 25), np.tile(stats["temp"]["min"], 25), np.tile(stats["disp"]["min"], 25), np.tile(stats["acceleration"]["min"], 25))), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.concatenate((np.tile(stats["voltage"]["max"], 25), np.tile(stats["temp"]["max"], 25), np.tile(stats["disp"]["max"], 25), np.tile(stats["acceleration"]["max"], 25))), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.concatenate((np.tile(stats["voltage"]["mean"], 25), np.tile(stats["temp"]["mean"], 25), np.tile(stats["disp"]["mean"], 25), np.tile(stats["acceleration"]["mean"], 25))), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.concatenate((np.tile(stats["voltage"]["std"], 25), np.tile(stats["temp"]["std"], 25), np.tile(stats["disp"]["std"], 25), np.tile(stats["acceleration"]["std"], 25))), device=self.dev, dtype=self.dtype) + + self.hw = self.calc_hw() + + def calc_hw(self): + m = 120 # mass of sensor network infrastructure e.g. concrete bridge, machine ... in tons + k = 10000. # stiffness, 10.000 as default value + # fs = 10 # [Hz] Sampling Frequency, according to Sampling Rate of the Sensors, to be adjusted + + # Mass matrix + M = np.eye(self.num_sensors) * m + _ndof = M.shape[0] # number of DOF (5): has to fit data shape + + # Stiffness matrix + K = np.zeros((_ndof, _ndof)) + for i in range(_ndof): + if i == 0: + K[i, i] = 1 + K[i, i + 1] = -1 + elif i == (_ndof - 1): + K[i, i] = 1 + K[i, i - 1] = -1 + else: + K[i, i] = 2 + K[i, i + 1] = -1 + K[i, i - 1] = -1 + + K *= k + + lam, FI = LA.eigh(K, b=M) # Solving eigen value problem + lam = np.abs(lam) + + fn = np.sqrt(lam) / (2 * np.pi) # Natural frequencies + + # Unity displacement normalised mode shapes + FI_1 = np.array([FI[:, k] / max(abs(FI[:, k])) for k in range(_ndof)]).T + # Ordering from smallest to largest + FI_1 = FI_1[:, np.argsort(fn)] + fn = np.sort(fn) + + M_M = FI_1.T @ M @ FI_1 # Modal mass + + xi = 0.02 # dampening ratio for all modes (2%) (for prestressed concrete), to be adjusted according to material + # Modal dampening + C_M = np.diag(np.array([2 * M_M[i, i] * xi * fn[i] * (2 * np.pi) for i in range(_ndof)])) + + C = LA.inv(FI_1.T) @ C_M @ LA.inv(FI_1) # Dampening matrix + + fx = M + C + K + + fw = np.fft.fft(fx) + + # Filter + # _b, _a = signal.butter(5, (49.9, 50.1), fs=fs, btype='bandpass') #net signal filtered out 0 < Wn < 1, doesn't work for ultralow frequency rates: 0 < (50.1 bzw. 49.9) / (fs * 0.5) < 1 sein (fs=10) + # filtdata = signal.filtfilt(_b, _a, data, axis=0) # filtered data + + hw = 1 / fw.T + return torch.tensor(hw, dtype=torch.complex64, device=self.dev) + + def forward(self, x, edge_index, weight, sim): + z = self.model(x, edge_index, weight, sim) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, acceleration, disp, cov = batch.x, batch.edge_index, batch.acc, batch.disp, batch.cov + if batch_idx == 0: + self.epoch_start_time = time.time() + + num_sensors = acceleration.shape[-1] + + # fix edge_index as it's not built with all nodes to save memory + e = np.vstack((np.repeat(np.arange(num_sensors), num_sensors), np.tile(np.arange(num_sensors), num_sensors))) + e = e[:, np.where(e[0] != e[1])[0]] + add = np.arange(0, batch.num_graphs) + add = np.repeat(add, e.shape[1]) + add *= num_sensors + e1 = np.tile(e, batch.num_graphs) + e2 = e1 + add + edge_index = torch.tensor(e2, device=self.dev, dtype=torch.int64) + + disp = disp.view((batch.num_graphs, 1, -1, disp.shape[-1])).repeat((1, num_sensors, 1, 1)) + + acceleration = acceleration.view((batch.num_graphs, 5, 5, acceleration.shape[-1])) + + cov = cov.view((acceleration.shape[0], acceleration.shape[1], cov.shape[-2], cov.shape[-1])) + sim = torch.mean(cov, dim=1) + sim = sim[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + + Sff = torch.fft.fft(cov.reshape((batch.num_graphs, -1, num_sensors, num_sensors)), axis=1) + + Syy = Sff @ self.hw.T + u, s, vh = torch.linalg.svd(Syy) + weights = u.mean(1) + + weights = weights[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + weights = torch.view_as_real(weights) + + # reshape and add acceleration to x + x = x.view((batch.num_graphs, 1, -1, x.shape[-1])) + x = x.repeat(1, num_sensors, 1, 1) + x = torch.concatenate((x.view((x.shape[0], x.shape[1], -1)), disp.view((x.shape[0], x.shape[1], -1)), + acceleration.permute(0, 3, 1, 2).reshape(acceleration.shape[0], acceleration.shape[3], + -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + # x = (x - self.mean) / self.std + x = (x - self.min) / (self.max - self.min) + + return x, edge_index, weights, sim + + def batch_cov(self, points): + """https://stackoverflow.com/a/71357620""" + B, N, D = points.size() + mean = points.mean(dim=1).unsqueeze(1) + diffs = (points - mean).reshape(B * N, D) + prods = torch.bmm(diffs.unsqueeze(2), diffs.unsqueeze(1)).reshape(B, N, D, D) + bcov = prods.sum(dim=1) / (N - 1) # Unbiased estimate + return bcov # (B, D, D) + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + losses = torch.nn.MSELoss()(y_hat, x.view(-1)) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/Luxemburg/LuxModeConvFastModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/LuxModeConvFastModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/LuxModeConvLaplaceModelPL.py b/ModeConvModel/models/Luxemburg/LuxModeConvLaplaceModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..afe65a29974d80744935477da1e285891651368c --- /dev/null +++ b/ModeConvModel/models/Luxemburg/LuxModeConvLaplaceModelPL.py @@ -0,0 +1,337 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy import linalg as LA +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.ModeConvLaplace import ModeConvLaplace +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics + + +class GNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvLaplace(14 * 25, args.bottleneck, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(14 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ModeConvLaplace(14 * 25, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(14 * 25, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, args.bottleneck, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index, weight, sim): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index, weight=weight, sim=sim) + lin(x) + x = x.relu() + return x + + +class EdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList( + [ModeConvLaplace(args.bottleneck, 14 * 25, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.conv_layers = nn.ModuleList( + [ModeConvLaplace(args.bottleneck, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append( + ModeConvLaplace(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, 14 * 25, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index, weight, sim): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index, weight=weight, sim=sim) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index, weight=weight, sim=sim) + self.lin_layers[-1](z) + return z.view(-1) + + +class GNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = GNNEncoder(num_sensors, args) + self.decoder = EdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index, weight, sim): + z = self.encoder(x, edge_index, weight, sim) + z = self.decoder(z, edge_index, weight, sim) + return z + + +class PSDLin(torch.nn.Module): + def __init__(self, num_sensors): + super().__init__() + self.lin_r = nn.Linear(num_sensors, num_sensors) + self.lin_i = nn.Linear(num_sensors, num_sensors) + + def forward(self, x): + z_r = self.lin_r(x) + z_i = self.lin_i(x) + return torch.complex(z_r, z_i) + + +class LuxModeConvLaplaceModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = GNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_lux/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.concatenate((np.tile(stats["voltage"]["min"], 25), np.tile(stats["temp"]["min"], 25), np.tile(stats["disp"]["min"], 25), np.tile(stats["acceleration"]["min"], 25))), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.concatenate((np.tile(stats["voltage"]["max"], 25), np.tile(stats["temp"]["max"], 25), np.tile(stats["disp"]["max"], 25), np.tile(stats["acceleration"]["max"], 25))), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.concatenate((np.tile(stats["voltage"]["mean"], 25), np.tile(stats["temp"]["mean"], 25), np.tile(stats["disp"]["mean"], 25), np.tile(stats["acceleration"]["mean"], 25))), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.concatenate((np.tile(stats["voltage"]["std"], 25), np.tile(stats["temp"]["std"], 25), np.tile(stats["disp"]["std"], 25), np.tile(stats["acceleration"]["std"], 25))), device=self.dev, dtype=self.dtype) + self.hw = self.calc_hw() + self.psd_lin = PSDLin(num_sensors) + + def calc_hw(self): + m = 120 # mass of sensor network infrastructure e.g. concrete bridge, machine ... in tons + k = 10000. # stiffness, 10.000 as default value + # fs = 10 # [Hz] Sampling Frequency, according to Sampling Rate of the Sensors, to be adjusted + + # Mass matrix + M = np.eye(self.num_sensors) * m + _ndof = M.shape[0] # number of DOF (5): has to fit data shape + + # Stiffness matrix + K = np.zeros((_ndof, _ndof)) + for i in range(_ndof): + if i == 0: + K[i, i] = 1 + K[i, i + 1] = -1 + elif i == (_ndof - 1): + K[i, i] = 1 + K[i, i - 1] = -1 + else: + K[i, i] = 2 + K[i, i + 1] = -1 + K[i, i - 1] = -1 + + K *= k + + lam, FI = LA.eigh(K, b=M) # Solving eigen value problem + lam = np.abs(lam) + + fn = np.sqrt(lam) / (2 * np.pi) # Natural frequencies + + # Unity displacement normalised mode shapes + FI_1 = np.array([FI[:, k] / max(abs(FI[:, k])) for k in range(_ndof)]).T + # Ordering from smallest to largest + FI_1 = FI_1[:, np.argsort(fn)] + fn = np.sort(fn) + + M_M = FI_1.T @ M @ FI_1 # Modal mass + + xi = 0.02 # dampening ratio for all modes (2%) (for prestressed concrete), to be adjusted according to material + # Modal dampening + C_M = np.diag(np.array([2 * M_M[i, i] * xi * fn[i] * (2 * np.pi) for i in range(_ndof)])) + + C = LA.inv(FI_1.T) @ C_M @ LA.inv(FI_1) # Dampening matrix + + fx = M + C + K + + fw = np.fft.fft(fx) + + # Filter + # _b, _a = signal.butter(5, (49.9, 50.1), fs=fs, btype='bandpass') #net signal filtered out 0 < Wn < 1, doesn't work for ultralow frequency rates: 0 < (50.1 bzw. 49.9) / (fs * 0.5) < 1 sein (fs=10) + # filtdata = signal.filtfilt(_b, _a, data, axis=0) # filtered data + + hw = 1 / fw.T + return torch.tensor(hw, dtype=torch.complex64, device=self.dev) + + def forward(self, x, edge_index, weight, sim): + # x = (x - self.mean) / self.std + z = self.model(x, edge_index, weight, sim) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def batch_cov(self, points): + """https://stackoverflow.com/a/71357620""" + B, N, D = points.size() + mean = points.mean(dim=1).unsqueeze(1) + diffs = (points - mean).reshape(B * N, D) + prods = torch.bmm(diffs.unsqueeze(2), diffs.unsqueeze(1)).reshape(B, N, D, D) + bcov = prods.sum(dim=1) / (N - 1) # Unbiased estimate + return bcov # (B, D, D) + + def process_batch(self, batch, batch_idx): + x, edge_index, acceleration, disp, cov = batch.x, batch.edge_index, batch.acc, batch.disp, batch.cov + if batch_idx == 0: + self.epoch_start_time = time.time() + + num_sensors = acceleration.shape[-1] + + e = np.vstack((np.repeat(np.arange(num_sensors), num_sensors), np.tile(np.arange(num_sensors), num_sensors))) + e = e[:, np.where(e[0] != e[1])[0]] + add = np.arange(0, batch.num_graphs) + add = np.repeat(add, e.shape[1]) + add *= num_sensors + e1 = np.tile(e, batch.num_graphs) + e2 = e1 + add + edge_index = torch.tensor(e2, device=self.dev, dtype=torch.int64) + + disp = disp.view((batch.num_graphs, 1, -1, disp.shape[-1])).repeat((1, num_sensors, 1, 1)) + + acceleration = acceleration.view((batch.num_graphs, 5, 5, acceleration.shape[-1])) + + cov = cov.view((acceleration.shape[0], acceleration.shape[1], cov.shape[-2], cov.shape[-1])) + sim = torch.mean(cov, dim=1) + sim = sim[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + + Sff = self.psd_lin(cov.reshape((batch.num_graphs, -1, num_sensors, num_sensors))) + + Syy = Sff @ self.hw.T + + weights = Syy.mean(1) + weights = weights[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + weights = torch.view_as_real(weights) + + # reshape and add acceleration to x + x = x.view((batch.num_graphs, 1, -1, x.shape[-1])) + x = x.repeat(1, num_sensors, 1, 1) + x = torch.concatenate((x.view((x.shape[0], x.shape[1], -1)), disp.view((x.shape[0], x.shape[1], -1)), + acceleration.permute(0, 3, 1, 2).reshape(acceleration.shape[0], acceleration.shape[3], + -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + # x = (x - self.mean) / self.std + x = (x - self.min) / (self.max - self.min) + + return x, edge_index, weights, sim + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + losses = torch.nn.MSELoss()(y_hat, x.view(-1)) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/Luxemburg/LuxModeConvLaplaceModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/LuxModeConvLaplaceModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/LuxMtGNNModelPL.py b/ModeConvModel/models/Luxemburg/LuxMtGNNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..fa6f33b502e48e8320c7f3a7f66c141bdc8c3cad --- /dev/null +++ b/ModeConvModel/models/Luxemburg/LuxMtGNNModelPL.py @@ -0,0 +1,197 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics +from Baseline_models.ssb_models.models.mtgnn.net import MtGNN + + +class MtDecoder(torch.nn.Module): + def __init__(self, args): + super().__init__() + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 14 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 14 * 25)) + + def forward(self, z): + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + return z + + +class MtGNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + self.encoder = MtGNN(gcn_depth=1, layers=args.num_layer, num_nodes=num_sensors, window_size=25, in_dim=14, out_dim=args.bottleneck, + skip_channels=args.hidden_dim, end_channels=args.hidden_dim, conv_channels=args.hidden_dim, residual_channels=args.hidden_dim) + if args.decoder == "linear": + self.decoder = MtDecoder(args) + elif args.decoder == "custom": + self.decoder = MtGNN(gcn_depth=1, layers=args.num_layer, num_nodes=num_sensors, window_size=1, + in_dim=args.bottleneck, out_dim=25 * 14, + skip_channels=args.hidden_dim, end_channels=args.hidden_dim, + conv_channels=args.hidden_dim, residual_channels=args.hidden_dim, subgraph_size=10) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, x): + z = self.encoder(x) + z = z.permute(0, 3, 2, 1) + z = self.decoder(z) + if self.args.decoder == "custom": + z = z.permute(0, 3, 2, 1) + return z + + +class LuxMtGNNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = MtGNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_lux/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.concatenate((stats["voltage"]["min"], stats["temp"]["min"], stats["disp"]["min"], [stats["acceleration"]["min"]])), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.concatenate((stats["voltage"]["max"], stats["temp"]["max"], stats["disp"]["max"], [stats["acceleration"]["max"]])), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.concatenate((stats["voltage"]["mean"], stats["temp"]["mean"], stats["disp"]["mean"], [stats["acceleration"]["mean"]])), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.concatenate((stats["voltage"]["std"], stats["temp"]["std"], stats["disp"]["std"], [stats["acceleration"]["std"]])), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, acceleration, disp = batch.x, batch.edge_index, batch.acc, batch.disp + if batch_idx == 0: + self.epoch_start_time = time.time() + acceleration = acceleration.view((batch.num_graphs, 5, 5, acceleration.shape[-1])).permute(0, 3, 1, 2).reshape( + (batch.num_graphs, 25, acceleration.shape[-1], 1)) + + num_sensors = acceleration.shape[-2] + disp = disp.view((batch.num_graphs, -1, 1, disp.shape[-1])).repeat((1, 1, num_sensors, 1)) + + x = x.view((batch.num_graphs, -1, 1, x.shape[-1])) + x = x.repeat(1, 1, num_sensors, 1) + x = torch.concatenate((x, disp, acceleration), dim=-1) + + # x = (x - self.mean) / self.std + x = (x - self.min) / (self.max - self.min) + + return x, edge_index + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index).view(train_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + x = x.permute(0, 2, 1, 3) + + losses = torch.nn.MSELoss()(y_hat, x) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y_hat = y_hat.view(val_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + + batchsize = len(val_batch.y) + x = x.permute(0, 2, 1, 3) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.reshape(-1).cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index).view(test_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + y = test_batch.y + + batchsize = len(test_batch.y) + x = x.permute(0, 2, 1, 3) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/Luxemburg/LuxMtGNNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/LuxMtGNNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/Luxemburg/__init__.py b/ModeConvModel/models/Luxemburg/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/ModeConvModel/models/Luxemburg/__init__.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/Luxemburg/__init__.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeAGCRNModelPL.py b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeAGCRNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..5e9a8d3099bc251c5b38c65c4e871cd08d0debdd --- /dev/null +++ b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeAGCRNModelPL.py @@ -0,0 +1,191 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics +from Baseline_models.ssb_models.models.agcrn.net import AGCRN + + +class AGCRNEdgeDecoder(torch.nn.Module): + def __init__(self, args): + super().__init__() + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + + def forward(self, z): + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + return z + + +class AGCRNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = AGCRN(num_sensors, input_dim=3, rnn_units=args.hidden_dim, output_dim=args.bottleneck, window_size=1, + dropout=0, num_layers=args.num_layer, embedding_dim=2, cheb_k=2) + if args.decoder == "linear": + self.decoder = AGCRNEdgeDecoder(args) + elif args.decoder == "custom": + self.decoder = AGCRN(num_sensors, input_dim=args.bottleneck, rnn_units=args.hidden_dim, output_dim=3*25, window_size=1, + dropout=0, num_layers=args.num_layer, embedding_dim=2, cheb_k=2) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, x): + z = self.encoder(x) + z = self.decoder(z) + return z + + +class SimulatedSmartBridgeAGCRNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = AGCRNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_simulated_smart_bridge/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.array([stats["inclination"]["min"], stats["temp"]["min"], stats["strain"]["min"]]), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.array([stats["inclination"]["max"], stats["temp"]["max"], stats["strain"]["max"]]), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.array([stats["inclination"]["mean"], stats["temp"]["mean"], stats["strain"]["mean"]]), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.array([stats["inclination"]["std"], stats["temp"]["std"], stats["strain"]["std"]]), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, strain = batch.x, batch.edge_index, batch.strain + if batch_idx == 0: + self.epoch_start_time = time.time() + + num_sensors = strain.shape[-1] + + strain = strain.view((batch.num_graphs, 5, 5, strain.shape[-1])) + + # reshape and add strain to x + x = torch.concatenate((x, strain.permute(0, 3, 1, 2).reshape(strain.shape[0] * strain.shape[3], -1)), dim=-1) + x = x.view((batch.num_graphs, num_sensors, 3, -1)).permute(0, 3, 1, 2) + + x = (x - self.min) / (self.max - self.min) + # x = (x - self.mean) / self.std + + return x, edge_index + + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + x = x.permute(0, 2, 1, 3).reshape(y_hat.shape) + + losses = torch.nn.MSELoss()(y_hat, x) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y_hat = y_hat.view(val_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + + batchsize = len(val_batch.y) + x = x.permute(0, 2, 1, 3) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.reshape(-1).cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y = test_batch.y + + batchsize = len(test_batch.y) + x = x.permute(0, 2, 1, 3).reshape(y_hat.shape) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeAGCRNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeAGCRNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeChebGNNModelPL.py b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeChebGNNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..1b3a50249ed04725bc976f0ba791df3bf80922ef --- /dev/null +++ b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeChebGNNModelPL.py @@ -0,0 +1,225 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn import ChebConv +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics + + +class ChebGNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ChebConv(3 * 25, args.bottleneck, K=5)]) + self.lin_layers = nn.ModuleList([Linear(3 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ChebConv(3 * 25, args.hidden_dim, K=5)]) + self.lin_layers = nn.ModuleList([Linear(3 * 25, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append(ChebConv(args.hidden_dim, args.hidden_dim, K=5)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ChebConv(args.hidden_dim, args.bottleneck, K=5)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index) + lin(x) + x = x.relu() + return x + + +class ChebEdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ChebConv(args.bottleneck, 3 * 25, K=5, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.conv_layers = nn.ModuleList( + [ChebConv(args.bottleneck, args.hidden_dim, K=5, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append( + ChebConv(args.hidden_dim, args.hidden_dim, K=5, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ChebConv(args.hidden_dim, 3 * 25, K=5, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index) + self.lin_layers[-1](z) + return z.view(-1) + + +class ChebGNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = ChebGNNEncoder(num_sensors, args) + self.decoder = ChebEdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index): + z = self.encoder(x, edge_index) + z = self.decoder(z, edge_index) + return z + + +class SimulatedSmartBridgeChebGNNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = ChebGNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_simulated_smart_bridge/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.repeat(np.array([stats["inclination"]["min"], stats["temp"]["min"], stats["strain"]["min"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.repeat(np.array([stats["inclination"]["max"], stats["temp"]["max"], stats["strain"]["max"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.repeat(np.array([stats["inclination"]["mean"], stats["temp"]["mean"], stats["strain"]["mean"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.repeat(np.array([stats["inclination"]["std"], stats["temp"]["std"], stats["strain"]["std"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x, edge_index) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, strain = batch.x, batch.edge_index, batch.strain + if batch_idx == 0: + self.epoch_start_time = time.time() + + strain = strain.view((batch.num_graphs, 5, 5, strain.shape[-1])) + + # reshape and add strain to x + x = torch.concatenate((x, strain.permute(0, 3, 1, 2).reshape(strain.shape[0] * strain.shape[3], -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + x = (x - self.min) / (self.max - self.min) + # x = (x - self.mean) / self.std + + return x, edge_index + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + + losses = torch.nn.MSELoss()(y_hat, x.view(-1)) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeChebGNNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeChebGNNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvFastModelPL.py b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvFastModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..bbe2723e1cb392714f4ba6df60f427ce887d9cdd --- /dev/null +++ b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvFastModelPL.py @@ -0,0 +1,323 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy import linalg as LA +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.ModeConvFast import ModeConvFast +from ModeConvModel.models.utils import compute_maha_threshold, compute_thresholds, generate_date_prefix, \ + compute_and_save_metrics + + +class GNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvFast(3 * 25, args.bottleneck, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(3 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ModeConvFast(3 * 25, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(3*25, args.hidden_dim)]) + for i in range(args.num_layer-2): + self.conv_layers.append(ModeConvFast(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvFast(args.hidden_dim, args.bottleneck, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index, weight, sim): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index, weight=weight, sim=sim) + lin(x) + x = x.relu() + return x + + +class EdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvFast(args.bottleneck, 3 * 25, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.conv_layers = nn.ModuleList([ModeConvFast(args.bottleneck, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer-2): + self.conv_layers.append(ModeConvFast(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvFast(args.hidden_dim, 3 * 25, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index, weight, sim): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index, weight=weight, sim=sim) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index, weight=weight, sim=sim) + self.lin_layers[-1](z) + return z.view(-1) + + +class GNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = GNNEncoder(num_sensors, args) + self.decoder = EdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index, weight, sim): + z = self.encoder(x, edge_index, weight, sim) + z = self.decoder(z, edge_index, weight, sim) + return z + + +class PSDLin(torch.nn.Module): + def __init__(self, num_sensors): + super().__init__() + self.lin_r = nn.Linear(num_sensors, num_sensors) + self.lin_i = nn.Linear(num_sensors, num_sensors) + + def forward(self, x): + z_r = self.lin_r(x) + z_i = self.lin_i(x) + return torch.complex(z_r, z_i) + + +class SimulatedSmartBridgeModeConvFastModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = GNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + self.cov_max = 370000 + + with open('processed_simulated_smart_bridge/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.repeat(np.array([stats["inclination"]["min"], stats["temp"]["min"], stats["strain"]["min"]]), [25,25,25]), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.repeat(np.array([stats["inclination"]["max"], stats["temp"]["max"], stats["strain"]["max"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.repeat(np.array([stats["inclination"]["mean"], stats["temp"]["mean"], stats["strain"]["mean"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.repeat(np.array([stats["inclination"]["std"], stats["temp"]["std"], stats["strain"]["std"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + + self.hw = self.calc_hw() + + def calc_hw(self): + m = 4377 # mass of sensor network infrastructure e.g. concrete bridge, machine ... in tons + k = 10000. # stiffness, 10.000 as default value + # fs = 10 # [Hz] Sampling Frequency, according to Sampling Rate of the Sensors, to be adjusted + + # Mass matrix + M = np.eye(self.num_sensors) * m + _ndof = M.shape[0] # number of DOF (5): has to fit data shape + + # Stiffness matrix + K = np.zeros((_ndof, _ndof)) + for i in range(_ndof): + if i == 0: + K[i, i] = 1 + K[i, i + 1] = -1 + elif i == (_ndof - 1): + K[i, i] = 1 + K[i, i - 1] = -1 + else: + K[i, i] = 2 + K[i, i + 1] = -1 + K[i, i - 1] = -1 + + K *= k + + lam, FI = LA.eigh(K, b=M) # Solving eigen value problem + lam = np.abs(lam) + + fn = np.sqrt(lam) / (2 * np.pi) # Natural frequencies + + # Unity displacement normalised mode shapes + FI_1 = np.array([FI[:, k] / max(abs(FI[:, k])) for k in range(_ndof)]).T + # Ordering from smallest to largest + FI_1 = FI_1[:, np.argsort(fn)] + fn = np.sort(fn) + + M_M = FI_1.T @ M @ FI_1 # Modal mass + + xi = 0.02 # dampening ratio for all modes (2%) (for prestressed concrete), to be adjusted according to material + # Modal dampening + C_M = np.diag(np.array([2 * M_M[i, i] * xi * fn[i] * (2 * np.pi) for i in range(_ndof)])) + + C = LA.inv(FI_1.T) @ C_M @ LA.inv(FI_1) # Dampening matrix + + fx = M + C + K + + fw = np.fft.fft(fx) + + # Filter + # _b, _a = signal.butter(5, (49.9, 50.1), fs=fs, btype='bandpass') #net signal filtered out 0 < Wn < 1, doesn't work for ultralow frequency rates: 0 < (50.1 bzw. 49.9) / (fs * 0.5) < 1 sein (fs=10) + # filtdata = signal.filtfilt(_b, _a, data, axis=0) # filtered data + + hw = 1 / fw.T + return torch.tensor(hw, dtype=torch.complex64, device=self.dev) + + def forward(self, x, edge_index, weight, sim): + z = self.model(x, edge_index, weight, sim) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def batch_cov(self, points): + """https://stackoverflow.com/a/71357620""" + B, N, D = points.size() + mean = points.mean(dim=1).unsqueeze(1) + diffs = (points - mean).reshape(B * N, D) + prods = torch.bmm(diffs.unsqueeze(2), diffs.unsqueeze(1)).reshape(B, N, D, D) + bcov = prods.sum(dim=1) / (N - 1) # Unbiased estimate + return bcov # (B, D, D) + + def process_batch(self, batch, batch_idx): + x, edge_index, strain, cov = batch.x, batch.edge_index, batch.strain, batch.cov + if batch_idx == 0: + self.epoch_start_time = time.time() + + num_sensors = strain.shape[-1] + + strain = strain.view((batch.num_graphs, 5, 5, strain.shape[-1])) + + cov = cov.view((strain.shape[0], strain.shape[1], cov.shape[-2], cov.shape[-1])) + sim = torch.mean(cov, dim=1) + sim = sim[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + + Sff = torch.fft.fft(cov.reshape((batch.num_graphs, -1, num_sensors, num_sensors)), axis=1) + + Syy = Sff @ self.hw.T + u, s, vh = torch.linalg.svd(Syy) + weights = u.mean(1) + + weights = weights[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + weights = torch.view_as_real(weights) + + # reshape and add strain to x + x = torch.concatenate((x, strain.permute(0, 3, 1, 2).reshape(strain.shape[0] * strain.shape[3], -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + x = (x - self.min) / (self.max - self.min) + # x = (x - self.mean) / self.std + + return x, edge_index, weights, sim + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)) + losses = losses.mean() + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvFastModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvFastModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvLaplaceModelPL.py b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvLaplaceModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..b72975d55d94f26788e21e043382ae91015117a3 --- /dev/null +++ b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvLaplaceModelPL.py @@ -0,0 +1,321 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy import linalg as LA +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.ModeConvLaplace import ModeConvLaplace +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics + + +class GNNEncoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvLaplace(3 * 25, args.bottleneck, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(3 * 25, args.bottleneck)]) + else: + self.conv_layers = nn.ModuleList([ModeConvLaplace(3 * 25, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(3*25, args.hidden_dim)]) + for i in range(args.num_layer-2): + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, args.bottleneck, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.bottleneck)) + + def forward(self, x, edge_index, weight, sim): + for conv, lin in zip(self.conv_layers, self.lin_layers): + x = conv(x, edge_index, weight=weight, sim=sim) + lin(x) + x = x.relu() + return x + + +class EdgeDecoder(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + if args.decoder == "linear": + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + elif args.decoder == "custom": + if args.num_layer <= 1: + self.conv_layers = nn.ModuleList([ModeConvLaplace(args.bottleneck, 3 * 25, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.conv_layers = nn.ModuleList( + [ModeConvLaplace(args.bottleneck, args.hidden_dim, K=1, num_sensors=num_sensors)]) + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.conv_layers.append( + ModeConvLaplace(args.hidden_dim, args.hidden_dim, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + + self.conv_layers.append(ModeConvLaplace(args.hidden_dim, 3 * 25, K=1, num_sensors=num_sensors)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, z, edge_index, weight, sim): + if self.args.decoder == "linear": + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + elif self.args.decoder == "custom": + for conv, lin in zip(self.conv_layers[:-1], self.lin_layers[:-1]): + z = conv(z, edge_index, weight=weight, sim=sim) + lin(z) + z = z.relu() + z = self.conv_layers[-1](z, edge_index, weight=weight, sim=sim) + self.lin_layers[-1](z) + return z.view(-1) + + +class GNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.encoder = GNNEncoder(num_sensors, args) + self.decoder = EdgeDecoder(num_sensors, args) + + def forward(self, x, edge_index, weight, sim): + z = self.encoder(x, edge_index, weight, sim) + z = self.decoder(z, edge_index, weight, sim) + return z + + +class PSDLin(torch.nn.Module): + def __init__(self, num_sensors): + super().__init__() + self.lin_r = nn.Linear(num_sensors, num_sensors) + self.lin_i = nn.Linear(num_sensors, num_sensors) + + def forward(self, x): + z_r = self.lin_r(x) + z_i = self.lin_i(x) + return torch.complex(z_r, z_i) + + +class SimulatedSmartBridgeModeConvLaplaceModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = GNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + self.cov_max = 370000 + with open('processed_simulated_smart_bridge/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.repeat(np.array([stats["inclination"]["min"], stats["temp"]["min"], stats["strain"]["min"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.repeat(np.array([stats["inclination"]["max"], stats["temp"]["max"], stats["strain"]["max"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.repeat(np.array([stats["inclination"]["mean"], stats["temp"]["mean"], stats["strain"]["mean"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.repeat(np.array([stats["inclination"]["std"], stats["temp"]["std"], stats["strain"]["std"]]), [25, 25, 25]), device=self.dev, dtype=self.dtype) + self.hw = self.calc_hw() + self.psd_lin = PSDLin(num_sensors) + + def calc_hw(self): + m = 4377 # mass of sensor network infrastructure e.g. concrete bridge, machine ... in tons + k = 10000. # stiffness, 10.000 as default value + # fs = 10 # [Hz] Sampling Frequency, according to Sampling Rate of the Sensors, to be adjusted + + # Mass matrix + M = np.eye(self.num_sensors) * m + _ndof = M.shape[0] # number of DOF (5): has to fit data shape + + # Stiffness matrix + K = np.zeros((_ndof, _ndof)) + for i in range(_ndof): + if i == 0: + K[i, i] = 1 + K[i, i + 1] = -1 + elif i == (_ndof - 1): + K[i, i] = 1 + K[i, i - 1] = -1 + else: + K[i, i] = 2 + K[i, i + 1] = -1 + K[i, i - 1] = -1 + + K *= k + + lam, FI = LA.eigh(K, b=M) # Solving eigen value problem + lam = np.abs(lam) + + fn = np.sqrt(lam) / (2 * np.pi) # Natural frequencies + + # Unity displacement normalised mode shapes + FI_1 = np.array([FI[:, k] / max(abs(FI[:, k])) for k in range(_ndof)]).T + # Ordering from smallest to largest + FI_1 = FI_1[:, np.argsort(fn)] + fn = np.sort(fn) + + M_M = FI_1.T @ M @ FI_1 # Modal mass + + xi = 0.02 # dampening ratio for all modes (2%) (for prestressed concrete), to be adjusted according to material + # Modal dampening + C_M = np.diag(np.array([2 * M_M[i, i] * xi * fn[i] * (2 * np.pi) for i in range(_ndof)])) + + C = LA.inv(FI_1.T) @ C_M @ LA.inv(FI_1) # Dampening matrix + + fx = M + C + K + + fw = np.fft.fft(fx) + + # Filter + # _b, _a = signal.butter(5, (49.9, 50.1), fs=fs, btype='bandpass') #net signal filtered out 0 < Wn < 1, doesn't work for ultralow frequency rates: 0 < (50.1 bzw. 49.9) / (fs * 0.5) < 1 sein (fs=10) + # filtdata = signal.filtfilt(_b, _a, data, axis=0) # filtered data + + hw = 1 / fw.T + return torch.tensor(hw, dtype=torch.complex64, device=self.dev) + + def forward(self, x, edge_index, weight, sim): + z = self.model(x, edge_index, weight, sim) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def batch_cov(self, points): + """https://stackoverflow.com/a/71357620""" + B, N, D = points.size() + mean = points.mean(dim=1).unsqueeze(1) + diffs = (points - mean).reshape(B * N, D) + prods = torch.bmm(diffs.unsqueeze(2), diffs.unsqueeze(1)).reshape(B, N, D, D) + bcov = prods.sum(dim=1) / (N - 1) # Unbiased estimate + return bcov # (B, D, D) + + def process_batch(self, batch, batch_idx): + x, edge_index, strain, cov = batch.x, batch.edge_index, batch.strain, batch.cov + if batch_idx == 0: + self.epoch_start_time = time.time() + num_sensors = strain.shape[-1] + + strain = strain.view((batch.num_graphs, 5, 5, strain.shape[-1])) + + cov = cov.view((strain.shape[0], strain.shape[1], cov.shape[-2], cov.shape[-1])) + sim = torch.mean(cov, dim=1) + sim = sim[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + + Sff = self.psd_lin(cov.reshape((batch.num_graphs, -1, num_sensors, num_sensors))) + + Syy = Sff @ self.hw.T + + weights = Syy.mean(1) + weights = weights[:, edge_index[0, :edge_index.shape[1] // batch.num_graphs], + edge_index[1, :edge_index.shape[1] // batch.num_graphs]] + weights = torch.view_as_real(weights) + + # reshape and add strain to x + x = torch.concatenate((x, strain.permute(0, 3, 1, 2).reshape(strain.shape[0] * strain.shape[3], -1)), dim=-1) + x = x.view((-1, x.shape[-1])) + + x = (x - self.min) / (self.max - self.min) + # x = (x - self.mean) / self.std + + return x, edge_index, weights, sim + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)) + losses = losses.mean() + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + + batchsize = len(val_batch.y) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index, weights, sim = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index, weights, sim) + y = test_batch.y + + batchsize = len(test_batch.y) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x.view(-1)).view(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvLaplaceModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeModeConvLaplaceModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeMtGNNModelPL.py b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeMtGNNModelPL.py new file mode 100644 index 0000000000000000000000000000000000000000..8ec34a090748412b985c5759480500160f66e100 --- /dev/null +++ b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeMtGNNModelPL.py @@ -0,0 +1,193 @@ +import json +import os +import random +import time +from pathlib import Path + +import numpy as np +import pytorch_lightning as pl +import torch +from scipy.spatial.distance import mahalanobis +from torch import nn +from torch.nn import Linear +from torch.optim import Adam +from torch_geometric.data import Data +from torch_geometric.nn.dense.linear import Linear + +from ModeConvModel.models.utils import compute_thresholds, compute_maha_threshold, generate_date_prefix, \ + compute_and_save_metrics +from Baseline_models.ssb_models.models.mtgnn.net import MtGNN + + +class MtDecoder(torch.nn.Module): + def __init__(self, args): + super().__init__() + if args.num_layer <= 1: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, 3 * 25)]) + else: + self.lin_layers = nn.ModuleList([Linear(args.bottleneck, args.hidden_dim)]) + for i in range(args.num_layer - 2): + self.lin_layers.append(Linear(args.hidden_dim, args.hidden_dim)) + self.lin_layers.append(Linear(args.hidden_dim, 3 * 25)) + + def forward(self, z): + for lin in self.lin_layers[:-1]: + z = lin(z).relu() + z = self.lin_layers[-1](z) + return z + + +class MtGNNModel(torch.nn.Module): + def __init__(self, num_sensors, args): + super().__init__() + self.args = args + self.encoder = MtGNN(gcn_depth=1, layers=args.num_layer, num_nodes=num_sensors, window_size=25, in_dim=3, out_dim=args.bottleneck, + skip_channels=args.hidden_dim, end_channels=args.hidden_dim, conv_channels=args.hidden_dim, residual_channels=args.hidden_dim, subgraph_size=10) + if args.decoder == "linear": + self.decoder = MtDecoder(args) + elif args.decoder == "custom": + self.decoder = MtGNN(gcn_depth=1, layers=args.num_layer, num_nodes=num_sensors, window_size=1, in_dim=args.bottleneck, out_dim=25*3, + skip_channels=args.hidden_dim, end_channels=args.hidden_dim, conv_channels=args.hidden_dim, residual_channels=args.hidden_dim, subgraph_size=10) + else: + raise NotImplementedError(f"Decoder {args.decoder} not implemented") + + def forward(self, x): + z = self.encoder(x) + z = z.permute(0, 3, 2, 1) + z = self.decoder(z) + if self.args.decoder == "custom": + z = z.permute(0, 3, 2, 1) + return z + + +class SimulatedSmartBridgeMtGNNModelPL(pl.LightningModule): + def __init__(self, num_sensors: int, args) -> None: + super().__init__() + self.model = MtGNNModel(num_sensors, args) + self.num_sensors = num_sensors + self.save_hyperparameters() + self.learning_rate = 1e-3 if args.lr == "auto" else float(args.lr) + self.batch_size = args.batch_size + prefix = generate_date_prefix() + random.seed(args.seed) + self.args = args + self.prefix = os.path.dirname(os.path.abspath(__file__)) + f'/../../results/{args.dataset}/{args.model}/{prefix}' + Path(self.prefix).mkdir(parents=True, exist_ok=True) + self.epoch_duration = 0 + self.val_losses = [] + self.val_logits = [] + self.val_labels = [] + self.test_data = [] + self.test_scores = [] + self.test_maha_scores = [] + self.test_logits = [] + self.test_labels = [] + self.dev = "cpu" if args.no_cuda else "cuda" + with open('processed_simulated_smart_bridge/processed/metadata.json', 'r') as fp: + stats = json.load(fp) + self.min = torch.tensor(np.array([stats["inclination"]["min"], stats["temp"]["min"], stats["strain"]["min"]]), device=self.dev, dtype=self.dtype) + self.max = torch.tensor(np.array([stats["inclination"]["max"], stats["temp"]["max"], stats["strain"]["max"]]), device=self.dev, dtype=self.dtype) + self.mean = torch.tensor(np.array([stats["inclination"]["mean"], stats["temp"]["mean"], stats["strain"]["mean"]]), device=self.dev, dtype=self.dtype) + self.std = torch.tensor(np.array([stats["inclination"]["std"], stats["temp"]["std"], stats["strain"]["std"]]), device=self.dev, dtype=self.dtype) + + def forward(self, x, edge_index): + z = self.model(x) + return z + + def configure_optimizers(self): + optimizer = Adam(self.parameters(), lr=self.learning_rate) + return optimizer + + def process_batch(self, batch, batch_idx): + x, edge_index, strain = batch.x, batch.edge_index, batch.strain + if batch_idx == 0: + self.epoch_start_time = time.time() + strain = strain.view((batch.num_graphs, 5, 5, strain.shape[-1])) + + num_sensors = strain.shape[-1] + + # reshape and add strain to x + x = torch.concatenate((x, strain.permute(0, 3, 1, 2).reshape(strain.shape[0] * strain.shape[3], -1)), dim=-1) + x = x.view((batch.num_graphs, num_sensors, 3, -1)).permute(0, 3, 1, 2) + + x = (x - self.min) / (self.max - self.min) + # x = (x - self.mean) / self.std + + return x, edge_index + + def training_step(self, train_batch: Data, batch_idx): + x, edge_index = self.process_batch(train_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + x = x.permute(0, 2, 1, 3).reshape(y_hat.shape) + + losses = torch.nn.MSELoss()(y_hat, x) + + return losses + + def training_epoch_end(self, outputs) -> None: + self.epoch_duration = time.time() - self.epoch_start_time + self.log("loss/train", float(np.mean([x["loss"].cpu() for x in outputs]))) + + def validation_step(self, val_batch: Data, batch_idx): + x, edge_index = self.process_batch(val_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y_hat = y_hat.view(val_batch.num_graphs, self.num_sensors, x.shape[1], x.shape[-1]) + + batchsize = len(val_batch.y) + x = x.permute(0, 2, 1, 3) + losses = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + self.val_losses.append(losses) + + y = val_batch.y + self.val_labels.append(y) + self.val_logits.append(y_hat.reshape(-1).cpu()) + return losses.mean() + + def validation_epoch_end(self, outputs) -> None: + self.val_losses = torch.hstack(self.val_losses) + val_losses = self.val_losses.cpu().detach().numpy() + self.val_labels = torch.hstack(self.val_labels) + val_labels = self.val_labels.cpu().detach().numpy() + self.val_logits = torch.hstack(self.val_logits) + val_logits = self.val_logits.cpu().detach().numpy() + + if not self.args.no_maha_threshold: + val_df, mean, cov = compute_maha_threshold(val_labels, val_logits) + + self.maha_thresh = val_df["Threshold"][val_df["bal_acc"].argmax()] + self.maha_mean = mean + self.maha_cov = cov + + val_df, max_val_loss = compute_thresholds(val_losses, val_labels) + self.max_val_loss = max_val_loss + self.best_threshold = val_df["Threshold"][val_df["bal_acc"].argmax()] + + self.val_labels = [] + self.val_losses = [] + self.val_logits = [] + self.log("val_f1", val_df["F1"][val_df["bal_acc"].argmax()]) + self.log("threshold", self.best_threshold) + self.log("max_score_threshold", self.max_val_loss) + return self.log("loss/val", float(np.mean([x.cpu() for x in outputs]))) + + def test_step(self, test_batch: Data, batch_idx): + x, edge_index = self.process_batch(test_batch, batch_idx) + + y_hat = self.forward(x, edge_index) + y = test_batch.y + + batchsize = len(test_batch.y) + x = x.permute(0, 2, 1, 3).reshape(y_hat.shape) + scores = torch.nn.MSELoss(reduction="none")(y_hat, x).reshape(batchsize, -1).mean(axis=1) + + self.test_scores.append(scores) + self.test_labels.append(y) + + if not self.args.no_maha_threshold: + maha_scores = np.array([mahalanobis(data, self.maha_mean, self.maha_cov) for data in y_hat.cpu().detach().numpy().reshape(test_batch.num_graphs, -1)]) + self.test_maha_scores.append(maha_scores) + + def test_epoch_end(self, outputs) -> None: + compute_and_save_metrics(self) \ No newline at end of file diff --git "a/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeMtGNNModelPL.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/SimulatedSmartBridgeMtGNNModelPL.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/SimulatedSmartBridge/__init__.py b/ModeConvModel/models/SimulatedSmartBridge/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/ModeConvModel/models/SimulatedSmartBridge/__init__.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/SimulatedSmartBridge/__init__.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/__init__.py b/ModeConvModel/models/__init__.py new file mode 100644 index 0000000000000000000000000000000000000000..6c1d80c4975d916b2898db4ae67425f9df8dab32 --- /dev/null +++ b/ModeConvModel/models/__init__.py @@ -0,0 +1,10 @@ +from ModeConvModel.models.Luxemburg.LuxAGCRNModelPL import LuxAGCRNModelPL +from ModeConvModel.models.Luxemburg.LuxChebGNNModelPL import LuxChebGNNModelPL +from ModeConvModel.models.Luxemburg.LuxModeConvFastModelPL import LuxModeConvFastModelPL +from ModeConvModel.models.Luxemburg.LuxModeConvLaplaceModelPL import LuxModeConvLaplaceModelPL +from ModeConvModel.models.Luxemburg.LuxMtGNNModelPL import LuxMtGNNModelPL +from ModeConvModel.models.SimulatedSmartBridge.SimulatedSmartBridgeModeConvFastModelPL import SimulatedSmartBridgeModeConvFastModelPL +from ModeConvModel.models.SimulatedSmartBridge.SimulatedSmartBridgeChebGNNModelPL import SimulatedSmartBridgeChebGNNModelPL +from ModeConvModel.models.SimulatedSmartBridge.SimulatedSmartBridgeMtGNNModelPL import SimulatedSmartBridgeMtGNNModelPL +from ModeConvModel.models.SimulatedSmartBridge.SimulatedSmartBridgeAGCRNModelPL import SimulatedSmartBridgeAGCRNModelPL +from ModeConvModel.models.SimulatedSmartBridge.SimulatedSmartBridgeModeConvLaplaceModelPL import SimulatedSmartBridgeModeConvLaplaceModelPL \ No newline at end of file diff --git "a/ModeConvModel/models/__init__.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/__init__.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/select_model.py b/ModeConvModel/models/select_model.py new file mode 100644 index 0000000000000000000000000000000000000000..7e6078a8307ca76612e12700465602bd52a96f96 --- /dev/null +++ b/ModeConvModel/models/select_model.py @@ -0,0 +1,40 @@ +from ModeConvModel.models import LuxModeConvFastModelPL, LuxModeConvLaplaceModelPL +from ModeConvModel.models import LuxChebGNNModelPL +from ModeConvModel.models import LuxAGCRNModelPL +from ModeConvModel.models import LuxMtGNNModelPL +from ModeConvModel.models import SimulatedSmartBridgeModeConvFastModelPL +from ModeConvModel.models import SimulatedSmartBridgeAGCRNModelPL +from ModeConvModel.models import SimulatedSmartBridgeMtGNNModelPL +from ModeConvModel.models import SimulatedSmartBridgeChebGNNModelPL +from ModeConvModel.models import SimulatedSmartBridgeModeConvLaplaceModelPL + + +def select_model(args): + if args.dataset == "luxemburg": + if args.model == "ModeConvFast": + return LuxModeConvFastModelPL(26, args) + elif args.model == "ModeConvLaplace": + return LuxModeConvLaplaceModelPL(26, args) + elif args.model == "ChebConv": + return LuxChebGNNModelPL(26, args) + elif args.model == "AGCRN": + return LuxAGCRNModelPL(26, args) + elif args.model == "MtGNN": + return LuxMtGNNModelPL(26, args) + else: + raise NotImplementedError(f"model {args.model} does not exist") + elif args.dataset == "simulated_smart_bridge": + if args.model == "ModeConvFast": + return SimulatedSmartBridgeModeConvFastModelPL(12, args) + elif args.model == "ModeConvLaplace": + return SimulatedSmartBridgeModeConvLaplaceModelPL(12, args) + elif args.model == "ChebConv": + return SimulatedSmartBridgeChebGNNModelPL(12, args) + elif args.model == "AGCRN": + return SimulatedSmartBridgeAGCRNModelPL(12, args) + elif args.model == "MtGNN": + return SimulatedSmartBridgeMtGNNModelPL(12, args) + else: + raise NotImplementedError(f"model {args.model} does not exist") + else: + raise NotImplementedError(f"dataset {args.dataset} does not exist") \ No newline at end of file diff --git "a/ModeConvModel/models/select_model.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/select_model.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/ModeConvModel/models/utils.py b/ModeConvModel/models/utils.py new file mode 100644 index 0000000000000000000000000000000000000000..d9ffdbc65900efdf20add0d8e8118e39fc516554 --- /dev/null +++ b/ModeConvModel/models/utils.py @@ -0,0 +1,185 @@ +import pickle +import random +import string + +import numpy as np +import pandas as pd +import torch +from scipy.spatial.distance import mahalanobis +import datetime +import time +from sklearn.metrics import confusion_matrix, balanced_accuracy_score, roc_auc_score +from scipy import linalg + + +def generate_date_prefix(random_letters=True) -> str: + out = f'{str(datetime.date.today())}_{datetime.datetime.now().hour}-{datetime.datetime.now().minute}' + if random_letters: + t = 1000 * time.time() + random.seed(int(t) % 2 ** 32) + out = f'{out}_{"".join(random.choices(string.ascii_uppercase, k=5))}' + return out + + +def compute_maha_threshold(val_labels, val_logits): + thresholds = [] + precisions = [] + recalls = [] + f1s = [] + bal_accs = [] + + ids = np.where(val_labels == 0) + + val_logits_normal = val_logits.reshape((len(val_labels), -1))[ids] + val_logits = val_logits.reshape((len(val_labels), -1)) + mean = val_logits_normal.mean(0) + cov = np.cov(val_logits_normal.T, val_logits_normal.T)[:val_logits_normal.shape[-1], :val_logits_normal.shape[-1]] + cov = linalg.pinv(cov) + maha_dist = np.array([mahalanobis(x, mean, cov) for x in val_logits_normal]) + all_maha_dist = np.array([mahalanobis(x, mean, cov) for x in val_logits]) + for i in range(10): + maha_thresh = np.quantile(maha_dist, 0.9 + i / 100) + predictions = np.where(all_maha_dist > maha_thresh, True, False) + + tn, fp, fn, tp = confusion_matrix(val_labels, predictions, labels=[0, 1]).ravel() + precision = tp / (tp + fp) + recall = tp / (tp + fn) + f1 = 2 * precision * recall / (precision + recall) + bal_acc = balanced_accuracy_score(val_labels, predictions) + + thresholds.append(maha_thresh) + precisions.append(precision) + recalls.append(recall) + f1s.append(f1) + bal_accs.append(bal_acc) + + val_df = pd.DataFrame(np.vstack((thresholds, precisions, recalls, f1s, bal_accs)).T, + columns=["Threshold", "Precision", "Recall", "F1", "bal_acc"]) + + return val_df, mean, cov + + +def compute_thresholds(val_losses, val_labels): + ids = np.where(val_labels == 0) + max_val_loss = val_losses[ids[0]].max() + + thresholds = [] + precisions = [] + recalls = [] + f1s = [] + bal_accs = [] + + for threshold in (np.arange(0, max_val_loss * 2, max_val_loss * 2 / 200)): + predictions = np.where(val_losses > threshold, True, False) + + tn, fp, fn, tp = confusion_matrix(val_labels, predictions, labels=[0, 1]).ravel() + precision = tp / (tp + fp) + recall = tp / (tp + fn) + f1 = 2 * precision * recall / (precision + recall) + bal_acc = balanced_accuracy_score(val_labels, predictions) + + thresholds.append(threshold) + precisions.append(precision) + recalls.append(recall) + f1s.append(f1) + bal_accs.append(bal_acc) + + thresholds = np.array(thresholds) + precisions = np.array(precisions) + recalls = np.array(recalls) + f1s = np.array(f1s) + bal_accs = np.array(bal_accs) + + val_df = pd.DataFrame(np.vstack((thresholds, precisions, recalls, f1s, bal_accs)).T, + columns=["Threshold", "Precision", "Recall", "F1", "bal_acc"]) + + return val_df, max_val_loss + + +def compute_and_save_metrics(model): + model.test_scores = torch.hstack(model.test_scores) + model.test_labels = torch.hstack(model.test_labels) + test_scores = model.test_scores.cpu().detach().numpy() + test_labels = model.test_labels.cpu().detach().numpy() + + predictions = np.where(test_scores > model.best_threshold, True, False) + + auc = roc_auc_score(test_labels, test_scores) + tn, fp, fn, tp = confusion_matrix(test_labels, predictions).ravel() + precision = tp / (tp + fp) + recall = tp / (tp + fn) + f1 = 2 * precision * recall / (precision + recall) + bal_acc = balanced_accuracy_score(test_labels, predictions) + + # predictions2 = np.where(test_scores > model.max_val_loss, True, False) + + # tn2, fp2, fn2, tp2 = confusion_matrix(test_labels, predictions2).ravel() + # precision2 = tp2 / (tp2 + fp2) + # recall2 = tp2 / (tp2 + fn2) + # f12 = 2 * precision2 * recall2 / (precision2 + recall2) + # bal_acc2 = balanced_accuracy_score(test_labels, predictions2) + + with open(model.prefix + '/metrics.txt', 'w') as f: + f.write(f"AUC: {auc}\n\n") + f.write(f"Best val balanced acc threshold:\n") + f.write(f"Precision: {precision}\n") + f.write(f"Recall: {recall}\n") + f.write(f"F1: {f1}\n") + f.write(f"bal_acc: {bal_acc}\n") + f.write(f"TP: {tp}\n") + f.write(f"TN: {tn}\n") + f.write(f"FP: {fp}\n") + f.write(f"FN: {fn}\n\n") + # f.write(f"Max normal val loss threshold:\n") + # f.write(f"Precision: {precision2}\n") + # f.write(f"Recall: {recall2}\n") + # f.write(f"F1: {f12}\n") + # f.write(f"bal_acc: {bal_acc2}\n") + # f.write(f"TP: {tp2}\n") + # f.write(f"TN: {tn2}\n") + # f.write(f"FP: {fp2}\n") + # f.write(f"FN: {fn2}\n\n") + + with open(model.prefix + '/scores.pkl', 'wb') as f: + pickle.dump(test_scores, f) + with open(model.prefix + '/labels.pkl', 'wb') as f: + pickle.dump(test_labels, f) + + model.log("metrics/test/prec", float(precision)) + model.log("metrics/test/recall", float(recall)) + model.log("metrics/test/f1", float(f1)) + model.log("metrics/test/bal_acc", float(bal_acc)) + model.log("metrics/test/auc", float(auc)) + + # model.log("metrics/test/prec2", float(precision2)) + # model.log("metrics/test/recall2", float(recall2)) + # model.log("metrics/test/f12", float(f12)) + # model.log("metrics/test/bal_acc2", float(bal_acc2)) + + if not model.args.no_maha_threshold: + test_maha_scores = np.hstack(model.test_maha_scores) + maha_predictions = np.where(test_maha_scores > model.maha_thresh, True, False) + auc_maha = roc_auc_score(test_labels, test_maha_scores) + tn_maha, fp_maha, fn_maha, tp_maha = confusion_matrix(test_labels, maha_predictions).ravel() + precision_maha = tp_maha / (tp_maha + fp_maha) + recall_maha = tp_maha / (tp_maha + fn_maha) + f1_maha = 2 * precision_maha * recall_maha / (precision_maha + recall_maha) + bal_acc_maha = balanced_accuracy_score(test_labels, maha_predictions) + + with open(model.prefix + '/metrics.txt', 'a') as f: + f.write(f"Mahalanobis threshold:\n") + f.write(f"AUC: {auc_maha}\n") + f.write(f"Precision: {precision_maha}\n") + f.write(f"Recall: {recall_maha}\n") + f.write(f"F1: {f1_maha}\n") + f.write(f"bal_acc: {bal_acc_maha}\n") + f.write(f"TP: {tp_maha}\n") + f.write(f"TN: {tn_maha}\n") + f.write(f"FP: {fp_maha}\n") + f.write(f"FN: {fn_maha}\n\n") + + model.log("metrics/test/prec_maha", float(precision_maha)) + model.log("metrics/test/recall_maha", float(recall_maha)) + model.log("metrics/test/f1_maha", float(f1_maha)) + model.log("metrics/test/bal_acc_maha", float(bal_acc_maha)) + model.log("metrics/test/auc_maha", float(auc_maha)) \ No newline at end of file diff --git "a/ModeConvModel/models/utils.py\357\200\272Zone.Identifier" "b/ModeConvModel/models/utils.py\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/data/SimulatedSmartBridge/simulated_smart_bridge_test.csv\357\200\272Zone.Identifier" "b/data/SimulatedSmartBridge/simulated_smart_bridge_test.csv\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git "a/data/SimulatedSmartBridge/simulated_smart_bridge_train.csv\357\200\272Zone.Identifier" "b/data/SimulatedSmartBridge/simulated_smart_bridge_train.csv\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/SimulatedSmartBridge/AGCRN/GNN.statedict b/data/pretrained/SimulatedSmartBridge/AGCRN/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..2d81ca4ec55b8cca8db6a78fcfd2979e40ec629a Binary files /dev/null and b/data/pretrained/SimulatedSmartBridge/AGCRN/GNN.statedict differ diff --git "a/data/pretrained/SimulatedSmartBridge/AGCRN/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/SimulatedSmartBridge/AGCRN/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/SimulatedSmartBridge/ChebConv/GNN.statedict b/data/pretrained/SimulatedSmartBridge/ChebConv/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..fe6554a6d814e132211ae5cf753918084818efbe Binary files /dev/null and b/data/pretrained/SimulatedSmartBridge/ChebConv/GNN.statedict differ diff --git "a/data/pretrained/SimulatedSmartBridge/ChebConv/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/SimulatedSmartBridge/ChebConv/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/SimulatedSmartBridge/ModeConvFast/GNN.statedict b/data/pretrained/SimulatedSmartBridge/ModeConvFast/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..6be48dbdf6336a6a8114ada5d0475e32b1dd623c Binary files /dev/null and b/data/pretrained/SimulatedSmartBridge/ModeConvFast/GNN.statedict differ diff --git "a/data/pretrained/SimulatedSmartBridge/ModeConvFast/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/SimulatedSmartBridge/ModeConvFast/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/SimulatedSmartBridge/ModeConvLaplace/GNN.statedict b/data/pretrained/SimulatedSmartBridge/ModeConvLaplace/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..d81ae5f8f2faf83ee42a8f410577a1f0f3054237 Binary files /dev/null and b/data/pretrained/SimulatedSmartBridge/ModeConvLaplace/GNN.statedict differ diff --git "a/data/pretrained/SimulatedSmartBridge/ModeConvLaplace/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/SimulatedSmartBridge/ModeConvLaplace/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/SimulatedSmartBridge/MtGNN/GNN.statedict b/data/pretrained/SimulatedSmartBridge/MtGNN/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..9be549709960cfe803a639579e46b9928faa7e2f Binary files /dev/null and b/data/pretrained/SimulatedSmartBridge/MtGNN/GNN.statedict differ diff --git "a/data/pretrained/SimulatedSmartBridge/MtGNN/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/SimulatedSmartBridge/MtGNN/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/luxemburg/AGCRN/GNN.statedict b/data/pretrained/luxemburg/AGCRN/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..47999bc9ee746ea9161375316eabb778ff921b90 Binary files /dev/null and b/data/pretrained/luxemburg/AGCRN/GNN.statedict differ diff --git "a/data/pretrained/luxemburg/AGCRN/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/luxemburg/AGCRN/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/luxemburg/ChebConv/GNN.statedict b/data/pretrained/luxemburg/ChebConv/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..5f6934418781497be929b59d9b99496cee858192 Binary files /dev/null and b/data/pretrained/luxemburg/ChebConv/GNN.statedict differ diff --git "a/data/pretrained/luxemburg/ChebConv/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/luxemburg/ChebConv/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/luxemburg/ModeConvFast/GNN.statedict b/data/pretrained/luxemburg/ModeConvFast/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..d81b5120d755a027d3de66abb842dfc4df675777 Binary files /dev/null and b/data/pretrained/luxemburg/ModeConvFast/GNN.statedict differ diff --git "a/data/pretrained/luxemburg/ModeConvFast/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/luxemburg/ModeConvFast/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/luxemburg/ModeConvLaplace/GNN.statedict b/data/pretrained/luxemburg/ModeConvLaplace/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..1a31fcab986feff75c564cfb8461751053e8f2aa Binary files /dev/null and b/data/pretrained/luxemburg/ModeConvLaplace/GNN.statedict differ diff --git "a/data/pretrained/luxemburg/ModeConvLaplace/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/luxemburg/ModeConvLaplace/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391 diff --git a/data/pretrained/luxemburg/MtGNN/GNN.statedict b/data/pretrained/luxemburg/MtGNN/GNN.statedict new file mode 100644 index 0000000000000000000000000000000000000000..cdcaa82e312a1943af94c92f577341b390358413 Binary files /dev/null and b/data/pretrained/luxemburg/MtGNN/GNN.statedict differ diff --git "a/data/pretrained/luxemburg/MtGNN/GNN.statedict\357\200\272Zone.Identifier" "b/data/pretrained/luxemburg/MtGNN/GNN.statedict\357\200\272Zone.Identifier" new file mode 100644 index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391