Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

cannot replicate the within-subject experimental results of BCI Competition IV dataset 2a using PyTorch #48

Open
wwwyz02 opened this issue Dec 15, 2023 · 0 comments

Comments

@wwwyz02
Copy link

wwwyz02 commented Dec 15, 2023

I cannot replicate the within-subject experimental results of BCI Competition IV dataset 2a using PyTorch
Here is my code for data processing

def getdata_A_T_22(filename,ids):
    raw_gdf = mne.io.read_raw_gdf(filename, stim_channel="auto", verbose='ERROR',
                                  exclude=(["EOG-left", "EOG-central", "EOG-right"]))
    #rename
    raw_gdf.rename_channels(
        {'EEG-Fz': 'Fz', 'EEG-0': 'FC3', 'EEG-1': 'FC1', 'EEG-2': 'FCz', 'EEG-3': 'FC2', 'EEG-4': 'FC4',
         'EEG-5': 'C5', 'EEG-C3': 'C3', 'EEG-6': 'C1', 'EEG-Cz': 'Cz', 'EEG-7': 'C2', 'EEG-C4': 'C4', 'EEG-8': 'C6',
         'EEG-9': 'CP3', 'EEG-10': 'CP1', 'EEG-11': 'CPz', 'EEG-12': 'CP2', 'EEG-13': 'CP4',
         'EEG-14': 'P1', 'EEG-15': 'Pz', 'EEG-16': 'P2', 'EEG-Pz': 'POz'})
    # Pre-load the data
    raw_gdf.load_data()
    #Bandpass filter 4-40Hz
    raw_gdf.filter(4.0, 40.0, fir_design="firwin", skip_by_annotation="edge")
    # Select data from 0.5s to 2.5s after the cue.  
    tmin, tmax =.5 , 2.5
    #  events_id
    events, events_id = mne.events_from_annotations(raw_gdf)
    event_id = dict({'769': 7, '770': 8, '771': 9, '772': 10})
    #get epoches
    epochs = mne.Epochs(raw_gdf, events, event_id, tmin, tmax, proj=True, baseline=None, preload=True)
    # print(epochs)
    labels = epochs.events[:, -1]
    data = epochs.get_data()
    # resample to 128Hz
    original_sampling_rate = 250
    target_sampling_rate = 128
    original_length = data.shape[-1]
    target_length = int(original_length * target_sampling_rate / original_sampling_rate)
    # resample    
    resampled_data = np.zeros((data.shape[0], data.shape[1], target_length))
    for i in range(data.shape[0]):
        for j in range(data.shape[1]):
            resampled_data[i, j, :] = resample(data[i, j, :], target_length)
    # print results
    print(f"orginal shape: {data.shape}, new shape: {resampled_data.shape}")
    data = resampled_data
    labels[labels == 7] = 0
    labels[labels == 8] = 1
    labels[labels == 9] = 2
    labels[labels == 10] = 3
    X = torch.from_numpy(data).unsqueeze(1)
    y = torch.from_numpy(labels).long()
    X, y = shuffle(X, y, random_state=42)  
    #Performing four-fold cross-validation
    num_folds = 4
    skf = StratifiedKFold(n_splits=num_folds, shuffle=True, random_state=42)
    folded_data = []
    for train_index, val_index in skf.split(X, y):
        X_train, X_val = X[train_index], X[val_index]
        y_train, y_val = y[train_index], y[val_index]
        #print(y_val)
        folded_data.append((X_train, y_train, X_val, y_val))
    # save
    for i, fold_data in enumerate(folded_data):
        X_train, y_train, X_val, y_val = fold_data
        print('shape',X_train.shape)
        print(f"Fold {i + 1}: Training Set Size = {len(X_train)}, Validation Set Size = {len(X_val)}")
        np.save(f'./2a_T_22s/train_data_a0{ids}_{i+1}.npy', X_train)
        np.save(f'./2a_T_22s/test_data_a0{ids}_{i+1}.npy', X_val)
        np.save(f'./2a_T_22s/train_labels_a0{ids}_{i+1}.npy', y_train)
        np.save(f'./2a_T_22s/test_labels_a0{ids}_{i+1}.npy', y_val)
for i in range(9):
    #The third subject's event differs
    if i==3:
        continue
    print(i+1)
    getdata_A_T_22(f"./BCICIV_2a_gdf/A0{i+1}T.gdf",i+1)

here is my code for train

class EEGNet(nn.Module):
    def __init__(self,Chans, T,clas, F1, D, F2 ):
        super(EEGNet, self).__init__()
        self.drop_out = 0.5
        self.block_1 = nn.Sequential(
            nn.Conv2d(
                in_channels=1,  
                out_channels=F1,  
                kernel_size=(1, 32),  
                bias=False,
                padding='same'
            ),  
            nn.BatchNorm2d(F1)  
        )
        self.block_2 = nn.Sequential(
            nn.Conv2d(
                in_channels=F1,  
                out_channels=D*F1,  
                kernel_size=(Chans, 1), 
                groups=F1,
                bias=False, 
                padding='valid'
            ),  
            nn.BatchNorm2d(D*F1),  
            nn.ELU(),
            nn.AvgPool2d((1, 4)),  
            nn.Dropout(self.drop_out)  
        )

        self.block_3 = nn.Sequential(
            nn.Conv2d(
                in_channels=D*F1,  
                out_channels=D*F1,  
                kernel_size=(1, 16),  
                groups=D*F1,
                bias=False,
                padding='same'
            ),  
            nn.Conv2d(
                in_channels=D*F1,  
                out_channels=F2,  
                kernel_size=(1, 1),  
                bias=False,
                padding='same'
            ),  
            nn.BatchNorm2d(F2), 
            nn.ELU(),
            nn.AvgPool2d((1, 8)),  
            nn.Dropout(self.drop_out)
        )

        self.out = nn.Linear((F2 * (T//32)), clas)

    def forward(self, x):
        x = self.block_1(x)
        x = self.block_2(x)
        x = self.block_3(x)

        x = x.view(x.size(0), -1)
        x = self.out(x)
        return x
device = torch.device("cuda")
creation=nn.CrossEntropyLoss()  
creation=creation.to(device)
learning_rate=3e-4
beta1=0.9
beta2=0.999
accs=np.zeros((9,4))
for sub in range(9):
    for fold in range(4):
        print(sub,fold)
        # load data
        loaded_train_data = np.load(f"./2a_T_22s/train_data_a0{sub+1}_{fold+1}.npy")
        #print(f"./2a_T_22s/train_data_a0{sub+1}_{fold+1}.npy")
        loaded_test_data = np.load(f'./2a_T_22s/test_data_a0{sub+1}_{fold+1}.npy')
        loaded_train_labels = np.load(f'./2a_T_22s/train_labels_a0{sub+1}_{fold+1}.npy')
        loaded_test_labels = np.load(f'./2a_T_22s/test_labels_a0{sub+1}_{fold+1}.npy')
        print(loaded_train_data.shape, loaded_train_labels.shape)
        #  PyTorch Tensor
        loaded_train_data_tensor = torch.Tensor(loaded_train_data)
        loaded_test_data_tensor = torch.Tensor(loaded_test_data)
        loaded_train_labels_tensor = torch.Tensor(loaded_train_labels).long()
        loaded_test_labels_tensor = torch.Tensor(loaded_test_labels).long()
        print(loaded_train_data_tensor.shape, loaded_train_labels_tensor.shape)
        # TensorDataset
        loaded_train_dataset = TensorDataset(loaded_train_data_tensor, loaded_train_labels_tensor)
        loaded_test_dataset = TensorDataset(loaded_test_data_tensor, loaded_test_labels_tensor)

        #  DataLoader
        batch_size = 32
        train_loader = DataLoader(loaded_train_dataset, batch_size=batch_size, shuffle=True)
        test_loader = DataLoader(loaded_test_dataset, batch_size=batch_size, shuffle=False)

        model = EEGNet(22,256,4,8, 2, 16)       
        optimizer=torch.optim.Adam(model.parameters(),lr=learning_rate)  
        model_save_path = './2a_T_22/models/'
        os.makedirs(model_save_path, exist_ok=True)
        EPOCHS = 300
        train_steps = 0
        model.cuda()
        test_steps = 0
        losses = []
        accuracies = []
        ftrain_label = torch.tensor(loaded_train_labels).to(device)

        for epoch in range(EPOCHS):  
            model.train() 
            epoch_train_loss = 0.0
            for train_batch, (train_image, train_label) in enumerate(train_loader):

                train_image, train_label = train_image.to(device), train_label.to(device)  
                train_predictions = model(train_image)

                batch_train_loss = creation(train_predictions, train_label)
                optimizer.zero_grad()  
                batch_train_loss.backward()
                optimizer.step()
                epoch_train_loss += batch_train_loss.item()
                train_steps += 1
                if train_steps % 500 == 0:
                    print("train{},train loss{}".format(train_steps, batch_train_loss.item()))
                if epoch%100==0:
                        model.eval()
                        with torch.no_grad():
                            predictions=model(loaded_test_data_tensor.to(device)).cpu()
                            predictions.numpy()
                            test_acc=(predictions.argmax(dim=1)==torch.tensor(loaded_test_labels )).sum()
                            print("test acc",test_acc/len(loaded_test_labels))
            losses.append(epoch_train_loss)
            train_predictions = model(loaded_train_data_tensor.to(device))

            accuracy = (torch.argmax(train_predictions, dim=1) == ftrain_label).float().mean().item()
            accuracies.append(accuracy)


        torch.save(model.state_dict(), model_save_path + f"modela0_{sub+1}{fold+1}.pth".format(epoch + 1))
        print("finish!")
        with torch.no_grad():
            predictions = model(loaded_test_data_tensor.to(device)).cpu()
            print(predictions[10])
            test_loss = creation(predictions, torch.tensor(loaded_test_labels))
            predictions.numpy()

            test_acc = (predictions.argmax(dim=1) == torch.tensor(loaded_test_labels)).sum()
            print("test acc", test_acc / len(loaded_test_labels))
            print("test loss{}".format(test_loss.item()))
            accs[sub][fold]=test_acc/len(loaded_test_labels)
average_per_experiment = np.mean(accs, axis=1)
# box
plt.boxplot(accs.T)

plt.title('Boxplot of Average per Experiment')
plt.xlabel('Experiment')
plt.ylabel('Average Value')
plt.show()

Here is my results
e9a08f98575f7f91dbd0848d74f399c

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant